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ABSTRACT 

This  thesis  involves  investigation  of  linear  filtering 
models  as  a  means  of  generating  texture  in  images.  Various 
autoregressive  filter  models  are  used  to  generate  various 
textures,  and  the  results  are  analyzed  to  determine  rela- 
tionships between  filter  parameters  and  texture  characteris- 
tics. A  two-dimensional  counterpart  to  the  autoregressive 
integrated  moving  average  (ARIMA)  model  from  one-dimensional 
time  series  analysis  theory  is  developed  and  tested  for 
texture  modeling  applications.  All  these  models  are  driven 
by  white  noise,  and  to  the  extent  that  real  images  can  be 
reproduced  this  way,  advantages  in  image  texture  transmis- 
sion could  be  realized.  Results  of  this  work  indicate  that 
the  purely  autoregressive  models  work  well  for  some  types  of 
image  textures,  but  that  for  the  textures  studied  the  ARIMA 
model  is  not  particularly  suitable. 
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I.   INTRODUCTION 

The  purpose  of  this  thesis  is  to  investigate  the  types 
and  quantity  of  image  textures  generated  using  a  two- 
dimensional  (2-D)  extension  of  the  Autoregressive  Integrated 
Moving  Average  (ARIMA)  model.  For  the  one-dimensional 
(e.g.,  time  series)  case,  the  theories  and  formulas  describ- 
ing this  model  are  outlined  in  Box  and  Jenkins  [Ref.  l:pp. 
85-103],  The  1-D  ARIMA  model  is  useful  when  the  time  series 
to  be  modeled  is  not  stationary  but  exhibits  some 
homogeneity  in  the  sense  that,  except  for  statistical 
differences  between  parts  of  the  time  series,  these 
different  parts  of  the  process  behave  similarly.  In  these 
cases  some  suitable  difference  of  the  process  may  be 
stationary,  and  hence  may  be  accurately  modeled  by  an 
Autoregressive  Moving  Average  (ARMA)  or  a  purely 
Autoregressive  (AR)  model.  The  resulting  stationary  time 
series  (generated  by  an  appropriate  ARMA  or  AR  filter  with 
white  noise  input)  is  applied  to  an  integra-tion  or 
summation  filter  (the  inverse  of  the  difference  operation) 
to  generate  the  original  nonstationary  time  series.  [Ref. 
l:p.  85]   Figure  1-1  shows  a  block  diagram  of  this  process. 

This  work  attempts  to  extend  these  concepts  to  two- 
dimensional  signal  processing.  In  order  to  simplify  the 
model,  the  moving   average   (MA)   portion   of   it   will   be 
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Block  Diagram  for  the  Autoregressive 
Integrated  Moving  Average  Model 


eliminated  (i.e.,  no  zeros  in  the  filter  Z  transform)  so 
that  only  purely  AR  models  will  be  considered  for  stationary 
image  generation.  The  procedures  for  modeling  image 
textures  using  AR  models  with  white  noise  input  are  well 
established  [Ref.  2:pp.  454-456],  However,  a  suitable  two- 
dimensional  difference  operation  and  its  inverse  must  be 
found  to  implement  the  concepts  outlined  above. 
The  research  is  divided  into  four  areas: 

1)  Investigation  of  the  various  types  of  image  textures 
generated  using  AR  models  where  filter  coefficients 
and  size  are  determined  a)  arbitrarily,  b)  using  a 
two-pole  separable  model,  and  c)  using  a  four-pole 
separable  model.  Separability  refers  to  the  fact  that 
the  Z  transform  of  the  AR  filter  can  be  factored  into 
components  representing  each  dimension  or  direction  of 
the  image. 

2)  Selection  of  a  difference  operator  and  a  realizable 
inverse  (integration  or  summation)  filter. 

3)  Application  of  the  above  autoregressively  generated 
images  to  the  summation  filter,  and  evaluation  of 
these  results. 

4)  Attempted  reproduction  of  actual  images  textures  using 
AR  models  whose  coefficients  are  determined  using  the 
statistics  of  the  image,  and  comparison  of  these 
results  to  those  obtained  by  a)  applying  the 
difference  operator  to  the  real  image,  b)  finding  the 
coefficients  of  the  AR  model  that  reproduce  the 
difference  image,  and  c)  applying  the  difference  image 
to  the  summation  filter.   This  comparison  was  intended 


to  discover  what  improvements,  if  any,  may  be  realized 
using  the  ARIMA  model  vice  a  purely  ARMA  (or  AR) 
model . 

The  remainder  of  the  thesis  is  organized  as  follows. 
Chapter  II  contains  methods  and  results  of  investigating 
various  autoregressive  image  models.  Chapter  III  deals  with 
the  formulation,  development,  and  testing  of  the  two- 
dimensional  summation  filter.  Chapter  IV  contains  the 
results  of  applying  various  AR-generated  images  to  the 
summation  filter  and  addresses  the  application  of  the  ARIMA 
model  to  real  image  textures.  Chapter  V  outlines  conclu- 
sions on  the  results  and  applicability  of  the  ARIMA  model. 
Although  the  ARIMA  modeling  was  not  highly  successful  ,  in 
reproducing  the  textures  studied  here,  plausible  reasons  are 
given  for  their  failure  and  conjectures  are  made  about  those 
circumstances  where  the  model  would  be  more  successful. 
Appendix  A  provides  information  on  the  computer  algorithms 
used  to  implement  the  equations  governing  the  above 
processes.  Appendices  B  through  G  contain  derivations  of 
spectral  and  autocorrelation  equations,  and  the  correspond- 
ing graphical  results,  governing  the  AR  processes  in  Chapter 
II.  Appendices  H  and  I  contain  graphical  results  associated 
with  the  inverse  filter  development  in  Chapter  III. 

Image  data  were  generated  using  computer  programs 
written  in  FORTRAN,  compiled  using  Version  4.5  under  the 
VAX/VMS  Version  4.4  operating  system.  The  images  were 
displayed  on  the  COMTAL  Vision  One/2  0.    The  gray  level 


intensity  range  of  pixel  values  is  0  (darkest)  to  255 
(lightest) ,  so  the  image  data  generated  had  to  be  scaled  to 
that  range  for  display  (see  Appendix  A) . 
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II.   THE  AUTOREGRESSIVE  IMAGE  MODEL 

A  two-dimensional  signal  (such  as  an  image  texture)  can 
be  modeled  using  a  two-dimensional  AR  model  with  white  noise 
input.  The  governing  equations  in  the  spatial  domain  are  of 
the  following  form  [Ref.  3:pp.  325-326]: 

P-1  Q-1 
Y(n,m)  =      -     I         I       aijy(n-i,m-j)  +  w(n,m)        (2.1) 
i=0  j=0 
(i,j)?^(0,0) 

where  y(n,m)  is  a  signal  representing  the  generated  image 
texture  at  pixel  location  (n,m)  ,  a-j^-;  is  the  filter  coeffi- 
cient matrix,  and  w(n,m)  is  a  two-dimensional  white  noise 
signal.  The  system  function  corresponding  to  the  filter  of 
Eq.  (2.1)  is  given  by 


1 


Y(Zi,Z2)   = =0 =1     -1  -1    7       ^- (P-1)  -(Q-1) 

1+^10^1  ■'^01^2  ^^11^1  ^2  "■•  •  '^^-iQ-lh  ^2 

(2.2) 
•  W(z-|_,Z2) 


where  z^  and  Z2  are  the  Z  transform  variables  corresponding 
to  spatial  coordinates  n  and  m.  Ideal  white  noise  has  an 
autocorrelation  function  that  is  an  impulse  and  a  flat  (con- 
stant) power  spectrum  with  magnitude  corresponding  to  the 
variance  of  the  white  noise  process  [Ref.  4:pp.  22-26]. 
Therefore  determination  of  the  filter  coefficients  a-j^j  will 
define  the  generated  image  process.    Procedures  will  be 
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outlined  later  to  estimate  the  coefficients  from  real  image 
data.  At  this  point  analytical  methods  will  be  used  to 
select  these  coefficients  and  the  resulting  images  will  be 
studied.  Figure  2-1  shows  an  example  of  a  white  noise  input 
image. 


Figure  2-1   White  Noise  Image 

A.   ARBITRARILY  SELECTED  FILTER  COEFFICIENTS;  P  =  2,  Q  =  2 

In  order  to  get  an  initial  idea  of  what  types  of  images 
might  be  generated  using  a  2  x  2  AR  filter  with  white  noise 
input,  filter  coefficients  were  at  first  selected  arbitrari- 
ly, but  subject  to  a  stability  constraint.  The  primary  con- 
straint on  coefficient  selection  is  that  of  filter 
stability.  Using  the  DeCarlo-Strintzis  Theorem  dealing  with 
multidimensional  filter  stability  [Ref.  3:pp.  197-198], 
alternately  setting  z^  =  1  and  Z2  =  1  and  determining  the 
location  of  the  pole  in  the  remaining  dimension  will 
indicate  whether  or  not  the  filter  is  stable.    If  the 
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magnitude  of  the  pole  in  the  remaining  dimension  is  less 
than  1,  the  filter  is  stable.  Even  with  this  condition, 
however,  there  are  an  infinite  number  of  possible  filter 
coefficient  combinations.  The  additional  constraint  of 
a]_o  ~  ^01  *-^^"  ^®  used,  and  comparisons  of  results  using 
various  values  of  a^^   can  be  made. 

Figure  2-2  shows  the  form  and  directionality  convention 
used  for  the  autoregressive  filter,  along  with  the  corres- 
ponding difference  equation  and  its  Z  transform. 


in(Z2) 


h(n,m) 


'01 


^11 

^10 


n(zi) 


1    1 
y(n,m)  =-     1       I     a.  . -y  (n-i,m-j)  +  w(n,m) 
i=0  j=0  ^^ 
(i,j)f^(0,0) 


Yiz^,z^)   =  H(z^,Z2) •W(z^,Z2)  = 


l+a^QZ-l4-aQ^Z2^+a^3_z^lz-^ 


W(z-|^,Z2) 


Figure  2-2   Autoregressive  Filter  Impulse  Response, 
Difference  Equation,  and  Z  Transform 


Although  it  is  difficult  to  make  precise  predictions  in 
two  dimensions,  one  can  expect  that  the  sign  and  magnitude 
of  aQ^o  o^  ^01  would  influence  the  correlation  between  pixels 
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in  the  corresponding  directions.  For  example,  a  positive 
value  for  a^^g  niight  be  expected  to  yield  an  image  with 
substantial  variation  in  the  n  direction  (low  correlation) , 
particularly  if  the  magnitude  of  a^g  i^  near  1.  A  negative 
value  for  agi  with  magnitude  near  1  might  yield  an  image 
with  lower  variations  in  pixel  intensity  (high  correlation) 
in  the  m  direction.  Since  the  filter  is  not  necessarily 
separable  (i.e.,  the  denominator  of  H( 2^,22)  cannot  be 
factored  into  the  form  D(Z2) • D(Z2) ) /  conclusions  drawn  from 
this  line  of  reasoning  may  not  be  completely  correct. 

Initial  attempts  at  generating  images  with  arbitrarily 
selected  coefficients  yielded  rather  uninteresting  results 
having  very  little  contrast  or  discernible  pattern. 
Continued  experimentation  with  combinations  where 
^10  ~  ^01  "^  ^  eventually  yielded  more  interesting  image  tex- 
tures. Figures  2-3  and  2-4  show  the  results  of  using  the 
constraint  a^o  =  ^oi  ~  -0.35  and  various  values  of  a^^i  for 
the  filter  coefficients.  For  positive  values  of  a^it  the 
images  are  rather  "grainy,"  with  higher  magnitudes  yielding 
a  somewhat  "finer"  graininess.  There  are  also  some  overall 
intensity  differences  observed.  For  the  negative  values  of 
^11'  the  results  are  much  more  interesting.  As  the  magni- 
tude increases,  there  is  a  gradually  more  noticeable  upper 
left  to  lower  right  orientation  of  the  image  texture,  and 
the  variations  from  lower  left  to  upper  right  become 
smoother   as   well.   Using  a^^Q  ~^1  ~  -0.38  and   a^i  =  -0.24 
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^10  ~  ^01  ~  -0.35 
a^i  =  0.8 


^10  ~  ^01  -    -0.35 
^11  ^  ^'^ 


^10 


=  am  =  -0.35 


'01 


^11 


=  0.2 


^10  ^  ^01  ^  -0.35 
^11  =0.0 


Sq^]^  —  SqI  ~  —0.35 


^11  ~  -0.1 


^10  ~  ^01  ~  -0.35 


^11  =  -0.2 


Figure  2-3   Images  Generated  Using  Arbitrarily- 
Selected  Filter  Coefficients 
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^10  ~  ^01  "=  -0.35 
^11  ~  -0.25 


^10  ~  ^01  =  -0.35 
^11  =  -0.27 


^10  ~  ^01  ~  -0.35 
^11  ~    "0-3 


^10  =  ^01  =  -0.38 
ail    ^  -0.24 


Figure  2-4   Images  Generated  Using  Arbitrarily 
Selected  Filter  Coefficients 
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yields  minor  variations  in  texture  pattern  and  overall  image 
intensity  when  compared  to  the  previous  case. 

Using  a^^g  ~  ^01  ~  0.35  and  varying  a^i  from  0.0  to  0.8 
(Figure  2-5) ,  the  images  obtained  deviate  very  little  from 
the  mean  intensity  value,  and  possess  minor  differences  in 
graininess.  With  these  positive  coefficients  some  negative 
correlation  might  be  expected,  and  the  fact  that  these 
images  are  "grainy"  indicates  the  existence  of  some  negative 
correlation  or  high  spatial  frequency  characteristics.  On 
initial  examination,  however,  the  low  contrast  of  the 
generated  images  tends  to  obscure  the  observed  graininess. 

B.   TWO  POLE,  SEPARABLE  AUTOREGRESSIVE  MODEL 

In  general,  it  is  difficult  to  relate  the  nature  or 
properties  of  a  two-dimensional  filter  to  the  precise  nature 
of  an  image  texture  that  may  be  generated  when  white  noise 
is  applied  to  that  filter.  In  order  to  simplify  the  effort 
and  to  obtain  a  better  understanding  of  the  problem,  the 
case  where  the  filter  (and  resulting  image  texture)  are 
separable  is  considered.  For  the  two  pole  separable  case 
considered  here,  the  filter  transfer  function  can  be 
factored  into  expressions  in  z-^  alone  and  Z2  alone.  The 
expressions  in  z^  and  Z2  each  have  one  pole  on  the  real  axis 
in  their  respective  Z  domains.  Figure  2-6  illustrates  the 
filter  structure,  the  corresponding  difference  equation,  and 
its  Z  transform. 
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^10  ~  ^01  ~  0.35 


^  \x 


<ifi 


^■.i^ 


5\  i**^*  *^*>*. 


^ 


^10  ~  ^01  =  0.35 
ail  =  0-2 


^10  ~  ^01  ~  0.35 
ail  ~  ^'^ 


^10  ~  ^01  ~  0.35 
ail  ~  ^'^ 


Figure  2-5   Images  Generated  Using  Arbitrarily 
Selected  Filter  Coefficients 
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h(n,in) 


ni(Z2) 


'01     ^11 
^10 


(ail  ■=  aio'aoi) 


n(zi) 


1    1 
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Figure  2-6   Autoregressive  Filter  Form,  Difference 
Equation,  and  Z  Transform 


Here  it  is  relatively  easy  to  relate  stability  of  the 
filter  to  the  location  of  the  poles  in  the  z^  and  Z2  planes. 
Since  the  quarter  plane  filter  is  separable  and  the 
components  are  causal,  one-dimensional  filter  stability 
theory  can  be  used  to  state  that  the  poles  in  each  plane 
must  have  magnitude  less  than  1  to  ensure  filter  stability. 
Figures  2-7  through  2-9  show  images  resulting  from  this 
model  for  various  values  of  a^Q  and  agi.  Note  that  the  sign 
convention  used  for  the  filter  difference  ecpiation  and  Z 
transform  results  in  poles  on  the  negative  side  of  the  real 
axis   for  positive  values  of  a-^Q   or  apj^,  and  vice  versa. 
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Figure  2-7   Images  Generated  Using  a  Two-Pole 
Autoregressive  Model 


20 


^10 


=  aoi  =  -0.25 


a^o  =  0.95   agi  =  0.5 


^10 


=0.5   aoi  =  0.95 


110  =  0-^5   aoi  =  0.75 


aio  =  -0.95   aoi  ^  0.25 


am  =  -0.95   am  =  -0.25 


^10 


01 


Figure  2-8   Images  Generated  Using  a  Two-Pole 
Autoregressive  Model 
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Figure  2-9   Images  Generated  Using  a  Two-Pole 
Autoregressive  Model 
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Careful  comparison  of  the  images  resulting  from  various  com- 
binations of  a^^o  ^^^  ^01  leads  to  the  following 
observations : 

1)  When  poles  are  located  in  the  same  place  in  the  z-^  and 
Z2  planes  on  the  negative  side  of  the  real  axis, 
magnitudes  near  1  yield  a  fine  graininess  with  patchy 
areas  and  low  overall  contrast.  As  the  magnitude  of 
the  pole  decreases,  the  graininess  becomes  more  coarse 
and  the  result  is  more  like  the  original  white  noise 
input.  No  directional  quality  in  the  image  pattern  is 
observed. 

2)  When  poles  are  located  in  the  same  place  on  the  posi- 
tive side  of  the  Z]^  and  Z2  real  axes,  a  somewhat  dif- 
ferent result  is  observed.  For  magnitudes  near  1,  an 
image  of  patchy  light  and  dark  areas  results,  with 
differing  amounts  of  correlation  between  pixels  in 
different  areas.  Slightly  discernible  "lines"  in  both 
the  horizontal  and  vertical  directions  are  also 
observed.  For  lower  pole  magnitudes  on  the  positive 
side  of  the  real  axis,  the  decreased  effect  of  the 
filter  on  the  white  noise  input  is  again  observed. 
This  result  is  more  like  the  white  noise  and  has  more 
contrast  than  the  corresponding  result  using  poles  on 
the  negative  side  of  the  real  axis. 

3)  For  pole  placements  in  the  Z]^  and  Z2  planes  which  are 
on  the  negative  side  of  the  real  axis  and  are  of 
unequal  magnitude,  the  results  have  a  very  fine 
graininess  and  low  contrast.  Some  slight  direction- 
ality is  observable  in  the  image  patterns,  with  lower 
frequency  variations  evident  in  the  direction  corres- 
ponding to  the  pole  with  smaller  magnitude. 

4)  For  pole  placements  in  the  Z]^  and  Z2  planes  which  are 
on  the  positive  side  of  the  real  axis  and  are  of 
unequal  magnitude,  much  more  directionality  and  varia- 
tion is  observable  in  the  image  pattern. 

5)  As  the  poles  are  placed  on  opposite  sides  of  the  real 
axis  and  are  separated  by  a  greater  distance,  direc- 
tionality becomes  more  evident  (with  higher  frequency 
variations  in  the  direction  of  the  more  negative 
pole) .  As  the  pole  separation  becomes  greater  and  as 
the  pole  magnitudes  become  closer  to  1,  smoother  sinu- 
soidal variations  are  evident. 

6)  When  the  pole  magnitudes  are  equal  and  have  opposite 
sign,  the  image  generated  using  pole  magnitudes  close 
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to  1  exhibits  high  frequency  sinusoidal  variations  in 
the  direction  of  the  negative  pole.  The  image  gener- 
ated with  the  lower  magnitude  poles,  as  would  be 
expected  from  the  above  results,  resembled  the 
unfiltered  white  noise. 

In  order  to  explain  these  image  patterns  analytically, 

analysis  of  the  power  spectrum  and  autocorrelation  function 

of  this  process  is  useful.   Since  this  model  is  separable, 

the  analysis  can  be  conducted  in  each  direction  separately. 

The  power  spectrum  is  defined  by  [Ref.  4:pp.  24-34]: 

Sv(w)   =  a2|H(eJ^)|2   =  a2H(e3^H(e-3^)  (2.3a) 


where 


H(e3W)   =   H(z)|2=e^^  (2.3b) 


and  for  this  case 


H(z)   =  ^3Y  (2.3c) 

1  +az 


Here  a^  is  the  magnitude  of  the  white  noise  power  spectrum. 
We  can  assume  that  a^  =  1  with  no  loss  of  generality  of  the 
results,  since  a^  does  not  affect  the  shape  of  the  frequency 
response. 

The  autocorrelation  function  is  related  to  the  filter 
transfer  function  through  the  equations  [Ref.  5: pp.  391- 
395]  : 
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h(n)   =   Z"l[H(z)]  (2.4a) 


y(n)   =   h(n)*w(n)  (2.4b) 


Ry(il)  =     o^       I       h(n)-h(n-£)  (2.4c) 


Specific  forms  of  the  power  spectrum  and  autocorrelation 
function  are  given  in  Appendix  B.  Since  Ry(£)  =  Ry(-2') 
[Ref.  5:p.  388],  calculating  the  expression  for  Ry(£),  ^-  <  0 
is  not  necessary.   From  Appendix  B,  the  results  are: 


Sy(aj)   =  2  (2.5) 

1  +2acos(w)  +a 

Ry(il)   =  -^^^     ^  >  0  (2.6) 

1  -a 


Appendix  C  shows  the  results  of  these  equations 
graphically  for  various  values  of  a  .  The  relationship  be- 
tween the  power  spectra  and  their  corresponding  autocorrela- 
tion functions  conforms  to  the  expected  results  from  theory 
(i.e.,  low  frequency  spectrum  with  smooth  autocorrelation 
function,  and  high  frequency  power  spectrum  with  rapidly 
varying  autocorrelation  function)  [Ref.  6:pp.  139-142].  The 
plots  in  Appendix  C  also  demonstrate  that:  1)  For  poles  on 
the  positive  side  of  the  real  axis  in  the  Z  plane  low  fre- 
quencies predominate  and  for  poles  on  the  negative  side  of 
the  real   axis  high   frequencies  predominate,  and  2)  Lower 
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magnitudes  of  a  result  in  a  broader  power  spectrum  and  a 
wider  range  of  frequencies  of  significant  magnitude.  Both 
of  these  observations  agree  with  the  image  results.  For 
images  generated  using  a  more  negative  pole  in  a  given 
direction,  fine,  high  frequency  graininess  is  observed  in 
that  direction  (though  the  low  contrast  or  low  variation 
about  the  mean  intensity  may  tend  to  make  this  effect  less 
noticeable) .  When  a  more  positive  pole  is  used,  lower  fre- 
quency variations  are  more  evident  in  the  corresponding 
direction.  As  lower  magnitudes  are  used  for  a  in  a  given 
direction,  more  random  variations  (indicative  of  a  wider 
range  of  significant  frequency  components)  are  observed  in 
that  direction.  The  form  of  the  autocorrelation  function 
for  these  cases  approaches  the  autocorrelation  function  for 
white  noise,  i.e.,  an  impulse.  Negative  poles  should  yield 
high  frequencies  since  the  negative  side  of  the  real  axis  in 
the  Z  plane  represents  a  digital  spatial  frequency  of  tt  , 
while  positive  poles  in  the  Z  plane  correspond  to  a  digital 
spatial  frequency  of  zero.  Note  that  even  when  the  poles 
are  placed  such  that  a  spatial  frequency  of  zero  should  pre- 
dominate, there  are  some  low  frequency  random  variations  in 
the  resulting  images.  Since  the  power  spectra  of  the  posi- 
tive poles  all  contain  some  non-zero  frequency  components 
(they  are  not  perfect  impulses  at  zero) ,  this  characteristic 
is  expected. 
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C.   FOUR  POLE,  SEPARABLE  AUTOREGRESSIVE  MODEL 

For  the  cases  considered  in  this  section,  H(z)  again  can 
be  factored  into  expressions  in  zjl  ^^id  Z2.  However  here 
each  factor  is  a  2nd  degree  polynomial  with  two  poles  in  the 
denominator.  Figure  2-10  illustrates  the  filter  structure, 
the  applicable  difference  equation,  and  the  corresponding  Z 
transform. 
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h(n,in) 


^02 
aoi 


^12 
^10 


^22 
321 
^20 


n(zi) 


2   2 

Y(n,m)  =  -     I       I       a.  .•y(n-i,m-j)  +w(n,m) 
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=1 =2  •  IT ^~^  •W(Zj_,Z2) 

1+a^QZ^  +^20^1    ^■'^01^2  -^^02^2 


-Vliz    .z   ) 

1^    -1^    -2^    -1^    -2^    -1  -1     -2  -1     -1  -2^    -2  -2  '^^''\'^2' 

^■'^lO^  -^^20^1  ■'^01^2^^02^2  ■'^ll^l  ^2  -^^21^1  ^2  ^^12^1  ^2  •'^22^1  ^2 


where 


^11  "  ^10*^01'   ^21  "  ^10*^02'   ^12  "  ^02 '^10'   ^22  ~  ^20*^02 


Figure  2-10   Autoregressive  Filter  Form,  Difference 
Equation,  and  Z  Transform 
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since  all  of  the  a-j^j  coefficients  are  real,  the  poles 
must  1)  both  be  on  the  real  axis,  or  2)  occur  in  complex 
conjugate  pairs  in  the  z-^   and  Z2  planes.   Again,  pole  magni- 
tudes must  be  less  than  1  to  ensure  filter  stability.   We 
will  assume  here  that  the  poles  in  each  of  the  factors  have 
equal  magnitudes  and  opposite  (or  0  or  ±  tt)  phase.   Letting 
a  1   =   magnitude  of  poles  in  the  z-^   plane 
9  2  =  pole  angle  (phase)  in  the  z-^   plane 
0.2   =   magnitude  of  poles  in  the  Z2  plane 
0  2=  pole  angle  (phase)  in  the  Z2  plane 
and  using  Euler's  relation,  the  denominators  of  H(Z]^)  and 
H(Z2)  can  be  expressed  as  follows: 


-1     -2         -'1  -1        -*  1  -1  -1  2  -2 

■'■"''^11^1   20^1  ^  (l-ci]_e   z^  )  (1-a^e    z^  )  =  l-2a^cos(0^) -z^  +a^z-|^ 


^^2  -1.   .  .     ,.  .  -1.  2  -2 


-1     -2        -'2  -1       ■'2   -1  -12- 

l+a^^z  +a  2^2  ^  ^'^-^2^       '^2   ^  ^^"^2^    ^2  ^  ^  l-2a200s(e2)  Z2  ■^^'2^2 


Hence: 


a^o  ~  -2oc-j^«cos(e  1)  aoi  =  -2a2*cos(62) 

2  _   2 

^20  =  °'l  ^02  ~  ^2 


This  gives  a  relationship  between  the  pole  magnitude  and 
angle  in  the  Z  domain  and  the  filter  coefficients  in  the 
spatial  domain. 
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1.   Complex  Conjugate  Poles 

Figures  2-11  through  2-13  illustrate  images 
generated  using  this  model  for  various  complex  conjugate 
pole  combinations  in  the  z^  and  Z2  directions.  In  comparing 
each  of  these  image  textures  in  terms  of  the  relative  effect 
of  pole  positioning  in  each  direction,  the  following  obser- 
vations can  be  made: 

1)  For  images  generated  using  poles  of  equal  magnitude 
and  angle  in  both  directions,  graininess  with  no 
directionality  to  the  pattern  resulted.  Higher  pole 
angles  yielded  finer  (higher  frequency)  graininess  and 
less  contrast.  Lower  magnitude  poles  yielded  a  more 
random  and  less  structured  graininess  pattern  at  the 
same  pole  angle. 

2)  Using  a  pole  angle  of  zero  (pole  on  positive  real 
axis)  in  one  direction  and  a  pole  of  some  non-zero 
angle  in  the  other  direction  yielded  images  with 
highly  directional  sinusoidal  patterns.  The  direction 
corresponding  to  the  pole  on  the  real  axis  was  not 
totally  devoid  of  variation,  but  variations  were  slow, 
i.e.,  of  very  low  frequency.  The  spatial  frequency  of 
the  sinusoidal  pattern  can  be  increased  by  increasing 
the  pole  angle.  Large  magnitude,  high  pole  angle 
combinations  yielded  much  cleaner  and  more  structured 
textures  than  low  magnitude,  low  pole  angle  combina- 
tions. Lower  magnitude,  high  pole  angle  combinations 
yielded  less  structured  textures  where  directionality 
was  evident  but  the  sinusoidal  pattern  was  obscured. 
Low  magnitude,  low  pole  angle  combinations  yielded 
very  random,  non-structured  textures  of  relatively 
high  contrast. 

3)  Using  poles  in  the  z-^  and  Z2  planes  with  the  same 
magnitude  but  different  pole  angles  resulted  in  some 
directionality  if  there  was  a  sufficiently  large 
magnitude  and  difference  in  the  pole  angles.  As 
observed  earlier,  the  direction  with  the  higher  pole 
angle  had  the  higher  spatial  frequency.  Large  pole 
magnitudes  (close  to  1)  resulted  in  more  structured 
but  rather  low-contrast  images  (the  low  contrast 
seemed  to  obscure  the  high  frequency  nature  of  the 
pattern  somewhat)  .  Pole  angles  in  Z]^  and  Z2  that  were 
of  close  value  made  it  difficult  to  detect  the  higher 
frequency    (higher   pole    angle)    direction.    Low 
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Figure  2-11  Images  Generated  Using  a  Four  Pole 
Autoregressive  Model  (Poles  Listed 
Below  Image) 
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Figure  2-12  Images  Generated  Using  a  Four  Pole 
Autoregressive  Model  (Poles  Listed 
Below  Image). 
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magnitudes  basically  negated  the  pole  angle  effects 
and  yielded  a  very  random,  unstructured,  high  contrast 
texture. 

4)  Where  pole  magnitudes  were  close  and  pole  angles  were 
different  in  the  z^  and  Z2  planes,  some  directionality 
in  the  texture  was  observed.  Again,  the  high  pole 
angle  direction  yielded  the  highest  frequency.  When 
the  pole  angles  in  both  Z  domains  had  similar  values 
and  the  magnitudes  of  the  poles  differed,  graininess 
with  little  or  no  discernible  directionality  or 
structure  resulted. 

These  observations  are  consistent  with  classical 
pole-zero  frequency  response  analysis  [Ref.  7:pp.  323-331]. 
There  is  a  direct  relationship  between  pole  angle  and 
spatial  frequency  in  a  given  direction,  and  the  magnitude  of 
the  poles  affects  the  amount  of  structure  and  definition  of 
the  sinusoidal  pattern  of  a  given  frequency  in  a  given 
direction.  Higher  magnitude  poles  result  in  a  narrower 
bandwidth  of  the  filter  and  yield  more  structure  and 
sinusoidal  pattern  definition.  Low  pole  magnitudes  give  the 
filter  wider  bandwidth  and  yield  images  with  less  structure 
and  definition  and  more  randomness  in  a  given  direction. 
While  directional  dependencies  are  evident  given  pole 
magnitude  and  angle  in  a  given  direction,  it  does  not  appear 
that  a  pattern  in  one  direction  is  totally  independent  of  a 
pattern  in  the  other  direction.  This  would  be  expected, 
even  though  the  model  is  separable,  due  to  the  cross  terms 
in  the  filter  structure. 

Filter  power  spectrum  and  autocorrelation  analysis 
can  be  conducted  in  this  case,  as  in  the  case  of  the  two 
pole  model.   The  derivations   for  Sy(w)   and  Ry(Jl)    are 
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somewhat  more  involved,  and  are  given  in  Appendix  D.   The 
resulting  expressions  for  Sy(aj)  and  Ry(^)  from  Appendix  D 


are: 


1 


^        l+a^+2(a^-2[a-^-a]oos(0)cos(co)+a  (cos  (2e)-K:os(2w) ) ) 


a^  ,cos{iQ)       cos((2+£)9)-a  cos(g,9).  ,.       q. 

Rv(^    =  T~^ ? 4 2 '       ^^-^'  (2.8) 

^        2sin^ei-a       1+a  -2a  cos(2e) 


The  plots  of  these  functions  for  the  various  pole 
magnitudes  and  angles  used  are  given  in  Appendix  E.  The 
9=0  case  is  equivalent  to  having  2  poles  on  the  real  axis 
at  a  given  magnitude  in  the  Z  plane.  As  would  be  expected, 
the  power  spectrum  for  each  model  showed  higher  magnitudes 
at  digital  frequencies  close  to  the  pole  angle.  Higher  pole 
magnitudes  yielded  sharper,  more  well-defined  power  spectrum 
magnitudes  at  the  given  frequency,  and  lower  pole  magnitudes 
yielded  less  well  defined  more  spread-out  power  spectra. 
Low  pole  magnitudes  almost  completely  obliterated  evidence 
of  low  frequency  power  spectrum  components,  and  degraded  its 
definition  and  sharpness  at  higher  frequencies.  This 
corresponds  to  the  observed  results  in  the  image  textures 
generated.  The  autocorrelation  functions  also  reflect  the 
appropriate  relationship  to  the  power  spectra  as  outlined  in 
the  discussion   for  the  two  pole   case,   i.e.,    greater 
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variation  in  the  autocorrelation  function  indicates  greater 
variation  between  pixels  a  given  distance  apart,  which  in 
turn  implies  higher  spatial  frequencies. 
2 .   Two  Real  Poles 

Rather  than  using  complex  conjugate  pole  locations 
to  obtain  real  filter  coefficient  values,  two  poles  on  the 
real  axis  may  also  be  used  for  a  given  direction.  They  may 
be  placed  at  different  locations  on  the  real  axis,  or  they 
may  be  placed  together.  The  latter  situation  is  equivalent 
to  the  6  =  0  (or  0  =  ±7t  if  placed  on  the  negative  real  axis) 
case,  as  mentioned  above.  For  the  two  real  pole  case,  the 
relevant  equations  are: 


H(z)   = 


TT        _l    ~  -1       -2 

(1-a  z  )  (1-a  z  )     l-(a^+a^)z  +a^o^z 


l+aQ3_z"^aQ2z"^ 


where 


^01   =   -(aa  +  ctb) 

^02   =  ^  a'^'b 

For  these  experiments  a  transfer  function  of  the 
complex  conjugate  pole  form  was  used  for  the  z^  direction, 
with  ai  =  0.9  and  G^  =  0.  For  the  Z2  direction  a  transfer 
function  with  two  poles  on  the   real   axis  was   used.   The 
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image  textures  that  result  for  various  values  of  a  ^  and  ctj^ 
are  given  in  Figures  2-14  and  2-15.  Some  observations  can 
be  made  about  these  results: 

1)  With  poles  placed  at  the  same  value  on  the  Z2  real 
axis,  rather  unstructured,  low  frequency  variations 
are  observed  in  the  Z2  direction.  The  more  positive 
poles  result  in  very  slow  variation  in  the  image 
texture,  while  the  lower  magnitude  positive  poles  show 
more  variation  in  the  Z2  direction.  When  the  poles 
are  placed  in  the  same  location  on  the  negative  side 
of  the  real  axis,  a  low  contrast  image  with  some 
noticeable  high  frequency  variations  results. 

2)  As  the  poles  are  moved  farther  apart  on  the  Z2  real 
axis,  high  frequency  variations  with  increasing  struc- 
ture and  oscillatory  form  are  evidenced  in  the  Z2 
direction. 

3)  When  poles  with  equal  magnitude  and  opposite  sign  are 
used,  fairly  structured  high  frequency  variations  are 
evidenced  in  certain  areas  of  the  image,  while  low 
frequency  variations  are  evident  in  other  areas  in  the 
Z2  direction.  Higher  magnitude  poles  yield  more 
discernible,  structured  variations,  while  lower 
magnitude  poles  of  opposite  sign  yield  discernible  but 
non-oscillatory  high  frequency  variations  in  certain 
areas  of  the  image. 

Of  particular  interest  is  the  fact  that  two  poles 

placed  at  the  same  value  on  the  negative  real  axis  in  the  Z2 

plane  yielded  some  high  frequency  variations.   This  is  in 

keeping  with  the  fact  that  values  on  the  negative  real  axis 

correspond   to   a   pole   angle   (and   corresponding   digital 

frequency)  of  6  =  tt  .   The  presence  of  poles  on  the  negative 

side  of  the  real  axis  of  the  Z2  plane  seems  to  give  rise  to 

the  high  frequency  variations  with  gradually  more  structure 

and  oscillatory  appearance  as  the  pole  is  moved  to  the  left 

(more  negative) . 
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aa  =  0.9   ab  =  0-9 


aa  =  0.5  a]-,  =  0.5 


^a  = 


-0.9   a)3  =  -0.9 


a^  =  0.9  aj^  =  o.5 


a^  =  0.9   a^  =  0. 0 


aa  =  0.9  a)3  =  -0.2 


Figure  2-14   Images  Generated  Using  a  Four  Pole 

(Two  Real  Poles)  Autoregressive  Model 
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aa  =  0-8  ab  =  -0.9 


a^=0.9  ai3  =  -0.9 


a^  =  0.5   aj^  =  -0.5 


Figure  2-15   Images  Generated  Using  a  Four  Pole 

(Two  Real  Poles)  Autoregressive  Model 
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The  expressions  for  the  power  spectrum  and  autocor- 
relation function  for  the  random  process  produced  by  driving 
the  filter  of  Figure  2-10  with  white  noise  are  derived  in 
Appendix  F.  It  is  shown  there  that  the  power  spectral 
density  and  the  autocorrelation  function  are  given  by: 


Sy  (co )   —  ^ 2  2 

1-2  (a  +0, +a  a, +a  a,  )cos(to)+2a  a,  cos  (2a)) +a 

2+z       i+1    ,  a+i       2+1 

,    a      a   a, +a  a,     a, 

Ry(M   =  1    [^   _^_^_^_b ^_D   ^  (2.10) 

■^         /     \2  T   2      1-a  0,       T  'i 
(V%^   ^-«a        ^^       ^~% 


Plots  of  these  functions  for  the  various  values  used 
in  this  section  are  given  in  Appendix  G.  The  power  spectrum 
results  are  consistent  with  the  observed  image  spatial  fre- 
quency characteristics.  Both  low  and  high  frequency 
components  were  contained  in  some  of  the  power  spectra,  and 
were  manifested  in  the  corresponding  images  as  both  low  and 
high  frequency  variations  in  the  Z2  direction.  The  nature 
of  the  autocorrelation  functions  related  to  the  power 
spectra  that  contained  low  and  high  frequency  components  was 
interesting.  Autocorrelation  functions  with  much  variation 
but  all  positive  values,  rather  than  the   equal   magnitude 
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positive  and  negative  values  evidenced  in  earlier  results, 
seems  to  reflect  the  higher  level  of  correlation  related  to 
the  low  frequency  (smoother  variations)  aspect  of  the  image 
texture  variations. 

D.   IMAGE  TEXTURE  ROTATION  TRANSFORMATION 

If  an  image  signal  x(n]^,n2)  consists  of  a  rotated 
version  of  another  image  w(m]^,m2)  such  that  n^^  =  Im^  +  Jm2 
and  n2  =  Km]^  +  Lm2 ,  where  I,  J,  K,  and  L  are  integers  and 
IL-KJ  7^  0,   then  the   Z   transform  X( 2^,22)   is   given  by 

T    K   T    L 
W(Z3^,  Z2,  z^,  Z2)  [Ref.  3:p.  182].   A  45°  rotation  corresponds 

to  1  =  1,  K=l,  J  =  l,  L=-l.    If  we  use  the  four  pole 


separable  result  for  H(Z]^,Z2),  as  shown  in  Figure  2-10,  and 
apply  the  above  transformation  {z-^  ->-  z^Z2,'z.2  -^  '^\'^2  )'  ^® 
find  after  simplification: 


Hr ( z ]_ z 2 )   =  _i  _-! _2  2     _]_  -1       ^^       -3-1 

l+aQ3_z^  Z2  +aQ2Zi  ^2'^^Q'^l  ^2   +^10^01^1  ''^10^02^1  ^2 

+  ^20^1^^2^^^.0^01^1^^2^''^20^02^1^  ^^'^^^ 


Notice  that  this  transfer  function  is  not  separable  but 
consists  of  a  rotated  version  of  a  separable  filter.  Figure 
2-16  illustrates  the  support  of  the  denominator  polynomial 
for  this  filter.  It  has  the  form  of  a  non-symmetric  half- 
plane  infinite  impulse  response  (IIR)  filter,  so  it  is 
recursively  computable. 
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^20^01 
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^20^02 

0 

^10^01 
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^10^02 
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^01 
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0 

— 

^02 
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Figure  2-16   Rotation  Transformation  Filter  Form 

The  application  of  this  filter,  using  filter  coefficients  of 
the  four  pole  separable  filter  with  poles  at  z^^  ^  0.9e  , 
Z2  ^  0.9e-3^  in  the  original  separable  filter  yielded  the 
result  shown  in  Figure  2-17. 


Figure  2-17   Result  of  Rotation  Transformation 

E .   SUMMARY 

Autoregressive  models  can  produce  a  variety  of  image 
textures.  For  general  two-dimensional  models,  the  system 
functions  are  generally  not   factorable   and   singularities 
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occur  on  surfaces,  not  at  isolated  points.  For  these 
reasons  it  is  difficult  to  design  two-dimensional  filters 
for  images  and  predict  the  resulting  character  of  the 
images.  Indeed,  even  to  ensure  stability  of  the  filter  is 
not  trivial.  As  a  result  we  concentrated  here  on  separable 
forms,  which  by  their  nature  are  much  easier  to  analyze. 
Certain  types  of  texture  patterns  using  various  separable 
autoregressive  models  can  be  predicted  based  on  filter  pole 
placement  in  the  z^  and  Z2  planes.  Arbitrary  or  random 
selection  of  filter  coefficients  can  yield  interesting  but 
generally  unpredictable  results.  Obviously,  an  infinite 
number  of  variations  on  the  models  above  could  be  attempted. 
Ultimately,  the  anticipated  utility  of  the  textures 
generated  will  guide  the  process  of  model  and  parameter 
selection. 


42 


III.   IMPLEMENTATION  OF  AN  FIR  SUMMATION  FILTER 

IN  TWO  DIMENSIONS 


To  implement  the  2-D  ARIMA  model,  the  inverse  of  a 
filter  representing  a  suitable  difference  operator  is 
needed.  One  possible  2-D  difference  operator  is  the 
Laplacian,  which  has  the  impulse  response  shown  in  Figure  3- 
1  [Ref.  8:pp.  212-213]. 


m 

h (n,m) 

0 

-1              0 

-1 

H    -1 

0 

-1               0 

n 


Figure  3-1   Laplacian  Impulse  Response 

Its  implementation  involves  convolving  it  with  an  image 
and  is  represented  by  the  following  difference  equation: 


1    1 
y(n,m)  =11     bi-ix(n-i,m-j ) 
i=-l  j=-l 


(3.1) 


where 


x(n,m)  is  the  image  input  signal,  and 

is  the  filter 
(bij  =  h(i,j)) 


b-j_j     is  the  filter  coefficient  matrix 
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In  the  Z-transform  domain  this  can  be  written  as: 

Y(Zi,Z2)   =   H(Zi,Z2) •X(Zi,Z2) 

=   (4  -  z"l  -  Zi  -  z~l  -  Z2)-X(Zi,Z2) 

In  areas  of  an  image  where  adjacent  pixels  have  similar  gray 
levels  (low  frequency,  homogeneous  areas) ,  the  result  of 
this  operator  will  be  approximately  zero.  Where  significant 
or  sharp  differences  in  gray  levels  between  adjacent  pixels 
exist,  the  result  of  this  operation  will  be  farther  from 
zero.  Thus  the  Laplacian  difference  operator  is  sometimes 
used  as  an  "edge  detector." 

The  problem  addressed  in  this  chapter  is  constructing 
the  inverse  of  the  operator.  In  the  Z  domain,  the 
expression  for  the  inverse  would  be  [Ref.  4: p.  36]: 


H-l(zi,Z2)   = 


H<^l-^2'      4-z-l-z,-z-^-Z3 


which  has  an  expansion  as  an  infinite  series  of  positive  and 
negative  powers  of  z-^  and  Z2.  That  is,  considering  this 
expression  as  a  problem  in  long  division,  the  result  of  such 
division  would  be  an  expression  of  the  form: 


H-l(zi,Z2)   =   .  Z   .1   a   z^ 

j_  =  — c»   -]=  — oo 


-^z:^ 


where 
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^ii  ~  coefficient  values  of  z"f-z3,  resulting  from  the 
long  division 


Note  that  if 


oo       oo 


1 =     y     y   a..zriz-j 

-1  -1  .  ^   .  ^    in  1  2 

then  cross  multiplication  would  yield: 


oo        00 

(4  -  z"l  -  Zi  -  z~2  -    Z2) ' _    1        .1   aijz"^z~^  =  1    (3.3) 

j^  =  — oo   -i=  — oo       -' 

The  double   summation  expression   in  z-^      and  Z2  will  be 

truncated  and  considered  to  be  an  FIR  filter  with  finite 

support  and  coefficients  aj^j  .   This  approximates  the  desired 

inverse  filter.    In  particular,  we  will  use  the  following 

constraints: 

1)   Choose  the  limits  of  summation  to  be  equal  in  both 
directions ,    i.e.. 


i=-L  j=-L   ^^  ^   2 


This  results  in  a  "square"  region  of  support  for  the 
filter  (all  values  outside  assumed  zero) . 

2)  Force  the  values  for  the  filter  coefficients  to  be 
symmetric,  i.e.,  a-j^-;  =  a-i-i  =  ^i--]  ~  ^-i-i  ~  ^i  i 
=  a_jj^  =  aj_j^  =  a_j_j^. 

Using   these   constraints   and   implementing   the   cross 

multiplication  equation  (3.3)  will  result  in  an  expression 

in  z^   and  Z2,  with  each  combination  of  the  z-'-z'^  terms  having 

a  coefficient  whose  form  is  a  summation  of  terms  in  aj^j 
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where  the  coefficients  of  a-j^-;  are  either  4  or  -1.  The 
coefficient  of  the  z?z^  term  must  equal  1  and  the  coeffi- 
cient of  any  other  zr-'-Z2^  term  must  equal  zero  to  satisfy 
equation  (3.3)  . 

As  an  example,  let  L  =  1.    Equation  (3.3)  can  then  be 
expressed  as: 

11 

(4    -    z~l    -    Zi    -    z~l    -    Z2)*        I  I     ai-;z~niz33 

i=-l  j=-l  -"J  -L   ^ 

1       1       11 
=   1  +  0-Zi  +  0'Z2  +  0*Z]_Z2  +  ... 


Performing  the  double  summation  yields 


l-yl      4-      a       _   -'^I'yO      4.      a .^l,,  — 1 


I  I       aijz-izj:3      =      a_i_izlzl   +   a_iozf-Z2 


+   a_T  T  ztz 


11' 


i=-l    j=-l 


1     2 


+  aooZ°z°  +  ao-izOzl  +  a^.iz-lz^  +  aiQZf ^ 


+  aoiZ2^  +  aiiz-^lzjl 


Performing  the  cross  multiplication  would  yield  45  different 
terms  in  various  combinations  of  z7^Z2^ .  Combining  these 
terms  to  find  the  coefficient  expression  for  each  z7^zj^ 
term  soon  becomes  rather  tedious  and  impractical  for  even 
moderate  values  of  L.  An  alternative  way  to  proceed  is  to 
choose  a  z^z^  term  on  the  right  hand  side  of  the  equation, 
and  for  each  term  in  the  expression  4-z'^-^-Z2-Z2"''-Z2  ^  deter- 
mine what  values  of  i  and  j  are  required  so  that  when 
each  term  is  multiplied  by  aijZ~^z~^,  it  will  result  in   an 
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expression  in  the  chosen  z?z^  term  on  the  right  hand  side  of 
the  equation.  Choosing  z^z^  (=  1)  on  the  right  hand  side  of 
the  equation,  we  have: 


4*ai-jZ,^Z2^  =  cz9zQ  when  i  =  0  and  j  =0;  so  c  =  4aoo 

-z~-^*a-i -;zr^z~^  =  cz*-'z'-'  when  i  =  -1  and  j  =0;  so  c  =  -a_Tn 

1  ^J  1   2       1  2  -^^ 

-Zi  •  a-i -;  zr^^z"^  =  cz^z*-*  when  i  =  1  and  j  =  0;  so  c  =  -an  a 
X   xj  ]^   2       1  2 

-z~-'--a-i -;  z~-'-z~3  =  cz*^z^  when  i  =  0  and  j  =-1;  so  c  =  -anj 

2  ID  1   2       1  2  "^ 

-Z9'an*-;z7"'"Z~^  =  cz'-'z'-'  when  i  =  0  and  j  =  1;  so  c  =  -am 

z   xj  2.   2       1  2 


So  the  coefficient  for  z^z^  is  simply  a  summation  of  the  c 

12 

terms  obtained  above,  i.e., 


(4aoo"  ^10  "  ^10  -  ^0-1  "  ^Ol)2°zO 

This  entire  expression  must  equal  1  to  satisfy  (3.3),  and 
.0^0  _ 


since  zp^  =  1,  4aoo  "  a.^Q  ~  ^10  "  ^0-1  "  ^01  =  1  also. 

Using  the  sa 
hand  side  yields: 


Using  the  same  method  for  the  z-^z-'-  term  on  the  right 

12 


4*ai-^z    ^z~^    =  cz?-z^  when  i   = -1   and  j    =-1;    so  c  =   4a_T_T 
-■■J    1212  -^  -'--'- 

-z~-^*a-i -;z~-^z~^    =   cz-^z^   when  i   =    -2    and   j    =   -1;    so   c   =   -a_9_T 
1^^1212  '^    ^ 

-Zi*aj^jz~^z~3    =  cz^z-^  when  i   =   0   and   j    =-1;    so   c  =   -ag-i 

-z~-^*a-i -;  z'-'-z"^    =   cz-'-z-'-   when  i   =   -1    and   j    =   -2;    so   c   =   -a_T_o 
2         ID     1212  -^^ 

-Zo 'a-i  ^  z~-'-z~3    =  cz-'-z-^  when  i   =   -1   and  j    =   0;    so   c  =   -a_nn 
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1-.1 


The  resulting  term  in  zfz-k    is: 


1^1 


(4a_i_l  -  a_2-i  -  ao-i  -  a_i_2  -  a_io)z^z 

This  expression  must  equal  zero,  since  there  is  no  z-'-z-'-  term 

on  the  right  side  of  (3.3),  so  4a_]^_]^  -  a_2_i  -  ao-1  ~  ^-1-2 

-  a_io  =  0. 

This  procedure  can  be  extended  to  any  number  of  zT-^z^ 

terms.   When  this  is  done,  the  resulting  expressions  for  the 

coefficients  of  z~^z~^  can  be  formed  into  a  set  of  simultan- 

1   2 

eous  equations  in  order  to  solve  for  the  a-^j  coefficient 
values.  However,  due  to  the  symmetry  condition  imposed 
above,  some  of  the  equations  for  the  coefficients  of  the 
z~^z~^  terms  are  linearly  dependent.  For  values  of  i  and  j 
that  yield  unique  or  distinct  values  for  a-j^j ,  the  resulting 
z~^z~^  coefficient  expressions  are  linearly  independent. 
For  example,  the  coefficient  expression  for  the  z~-^z^  term 
is  linearly  independent  of  the  coefficient  expression  for 
the  z'-'-z"?  term,  since  a^^g  ^  ^11*  ^^"^  "^^®  coefficient  ex- 
pression for  the  z-'-z^  term  is  linearly  dependent  on  the 
coefficient  expression  for  the  z'-'z"-^  term,  since  a^g  ~  ^01* 
Using  only  the  linearly  independent  equations  for  a  given 
filter  size  yields  a  set  of  p  equations  in  p  unknowns,  where 
p  is  the  number  of  unique  and  distinct  aj^j  values  in  the  in- 
verse filter.  The  value  of  p  is  related  to  the  size  of  the 
desired  inverse  filter.   If  the  size  of  the  filter  is   N  x  n. 
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the  number  of  unique  aj^j  values  using  the  symmetry  of 
constraint  above  is: 


p   =   (Nzl+i)  +  (Nil)  +  (Nil  _i)  +  ...  +  1  (3.4^) 


(N+1)  (N4-3) 
Z (3.4b) 


For  example,  with  N  =  7,  the  unique  a-j,  j  values  in  a  7x7 
inverse  filter  can  be  represented  by  a33,  ^32,  ^31/  ^30' 
^22'  ^21'  ^20'  ^11'  ^10'  ^00'  Though  there  are  49  elements 
in  a  7x7  filter,  all  of  them  are  equal  to  one  of  these 
values  listed,  due  to  symmetry.  Obviously,  N  is  constrained 
to  be  odd  for  a  square  filter  with  a  unique  element  ago  i^i 
the  center. 

The  solution  of  the  resulting  p  equations  yields  the 
values  for  the  p  filter  elements  or  coefficients.  This 
defines  the  FIR  approximation  to  the  inverse  filter.  It  is 
only  an  approximation  due  to  the  finite  size  constraint 
imposed,  and  it  might  be  expected  that  the  larger  the  filter 
size,  the  better  the  approximation. 

An  algorithmic  procedure  for  obtaining  the  a-j^j 
coefficients  is  outlined  below.   An  example  follows. 

1)  Determine  the  desired  size  of  the  inverse  filter. 

2)  For  each  combination  of  (positive)  i,j  values  corres- 
ponding to  a  unique  a-j^j  filter  coefficient,  identify 
th?  f  j.ve  term  summation  equation  associated  with  each 
z"^^Z2^  term. 
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3)  Combine  equal  aj^j  values  and  develop  a  matrix  of 
coefficients  for  the  aj^j  values.  Let  this  matrix  be 
A. 

4)  Denoting  the  vector  of  unique  a-j^j  values  as  a,  the  set 
of  simultaneous  equations  in  matrix  form  is: 


Aa 


rr  ^ 


(3.5) 


where 


[^00   ^10   ^11   ^20   ^21   ^22  •••  ^-1  N-1 


2   2 

The  top  row  of  A  corresponds  to  the  summation  of  aj^-; 
terms  that  represents  the  coefficient  of  the  zO"z! 
term. 

5)   Solve  (3.5)  for  a. 

An  example  of  this  procedure  is  appropriate  at  this 

point.   For  an  inverse  filter  of  size  N  x  N: 

Step  1 

Let  N  =  5  (therefore  L  =  2) 

Step  2 

The  coefficients  corresponding  to  each  unique  z~-'-z~^ 

term  are: 


zOzO 
12 

z-lzO 
1   2 

z-lz-1 
1  2 


"i^OO  ~  ^10  "  ^-10  ~  ^01  "  ^0-1 
4aio  ~  ^20  ~  ^00  "  ^11  "  ^1-1 
4an  -  a2i  -  agi  -  ^^2    "  ^^q 
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^1  ^2 

-> 

4a20 

z-2z-l 
1       2 

-> 

4a2i 

z-2z-2 

-> 

4a22 

^30  "  ^10  "  ^21  ~  ^2-1 
^31  ~  ^11  "  ^20  ~  ^2-2 
^32  ~  ^12  ~  ^23  ~  ^21 


Step  3 

Combining  equal  terms  in  Step  2  and  expressing  the 
coefficients  of  aj^-;  in  matrix  form  yields  an  A  matrix  of: 


A  = 


4 

-4 

0 

0 

0 

0 

1 

4 

-2 

-1 

0 

0 

0 

-2 

4 

0 

-2 

0 

0 

-1 

0 

4 

-2 

0 

0 

0 

-1 

-1 

4 

-1 

0 

0 

0 

0 

-2 

4 

Step  4 
With  A  given  in  Step  3,  the  a  vector  for  (3.5)  is 

a   =   [^00   ^10   ^11   ^20   ^21   ^22^ 

Step  5 

Solving  (3.5)  for  a  involves  inverting  A  and  multiplying 
itby[l   0   0   0   0   0]T.   Thus: 


a   =   A"l*[l   0   0   0   0   0]T 
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Appendix  H  illustrates  the  forms  of  the  resulting 
inverse  filters  of  various  sizes,  as  well  as  the  normalized 
and  unnormalized  filter  cross  sections. 

One  way  to  validate  the  resulting  inverse  filter  is  to 
convolve  it  with  the  original  Laplacian  difference  operator. 
The  result  should  approximate  an  impulse  at  the  origin. 
Appendix  I  shows  the  results  of  this  convolution  using  3x3, 
5x5,  7x7^  9x9,  15  X  15^  and  21  x  21  size  inverse  filters. 
It  is  seen  there  that  as  the  size  of  the  filter  gets  larger, 
it  becomes  a  better  approximation  to  the  true  inverse  and 
the  convolution  looks  more  like  an  impulse. 

To  test  the  application  of  this  filter  on  an  actual 
image,  a  test  image  was  filtered  using  the  Laplacian 
difference  operator.  Then  the  resulting  image  data  were 
filtered  again  using  various  size  inverse  filters.  The 
results  are  shown  in  Figures  3-2  and  3-3.  Note  that  the 
image  resulting  from  Laplacian  FIR  filtering  seems  more 
stationary  than  the  test  image,  which  was  one  of  the  desired 
results.  Inverse  filtering  of  that  result  yields  images 
that  are  progressively  more  similar  to  the  original  test 
image  as  the  size  of  the  inverse  filter  increases.  However, 
a  rather  large  inverse  filter  is  needed  to  accurately 
reproduce  the  image.  The  result  of  the  21  x  21  inverse 
filter  is  quite  similar  to  the  test  image,  with  some  loss  of 
contrast  or  darkness  in  certain  areas,  but  with  essentially 
the  same  pattern.   The  effect  of  the  size  limitation  of   the 
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Test  image 


Laplacian  filtered  test 
image 


3x3  inverse  filter 


5x5  inverse  filter 


Figure  3-2   Results  of  Filtering  Test  Image  with 
Laplacian  FIR  Filter  and  its  Inverse 
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7x7  inverse  filter 


9x9  inverse  filter 


15  X 15  inverse  filter 


21  X  21  inverse  filter 


Figure  3-3   Results  of  Filtering  Test  Image  with 
Laplacian  FIR  Filter  and  its  Inverse 
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inverse  filter,  as  manifested  in  the  convolution  of  the 
Laplacian  and  its  inverse  above,  would  seem  to  explain  the 
lack  of  perfect  test  image  reproduction.  Larger  inverse 
filter  sizes  could  be  tried,  but  large  inverse  filter  sizes 
relative  to  image  size  would  result  in  a  significant  portion 
around  the  edge  of  the  image  having  only  a  part  of  the 
filter  applied  to  it.  This  would  adversely  affect  the 
overall  quality  of  image  reproduction. 
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IV.   APPLICATION  OF  THE  ARIMA  MODEL  TO  IMAGE  TEXTURES 

As  outlined  in  Chapter  I,  the  utility  of  the  ARIMA  model 
centers  around  the  fact  that  a  difference  operator  applied 
to  an  image  texture  may  improve  the  stationarity  of  the 
image  statistical  characteristics.  A  stationary  image 
texture  is  required  for  accurate  modeling  by  autoregressive 
techniques,  and  it  was  hoped  that  application  of  the 
autoregressively  generated  texture  to  an  approximate  inverse 
of  the  difference  operator  may  yield  a  more  accurate  or 
recognizable  representation  of  the  original  nonstationary 
image,  as  compared  to  a  purely  autoregressively  generated 
version. 

A.   APPLICATION  OF  LAPLACIAN  INVERSE  FILTER  TO  AUTO- 
REGRESSIVELY GENERATED  IMAGES 

As  an  initial  examination  of  the  effects  of  the  inverse 

filter  developed  in  Chapter  III  on  image  textures,  selected 

images  generated  in  Chapter  II  were  input  to  the  21X21 

version  of  that  filter.   Figures  4-1  through  4-3  illustrate 

the  results.   All  attempts  resulted  in  a  blurred  or  smoothed 

version  of  the  original  image.    Since  the  inverse  of  a 

difference   operation   is   a   summation   or   "integration" 

operation,  and  since  integration  operations  can  be  expected 

to  smooth  or  blur  (low  pass  filter)  signals  [Ref.  8:pp.  136- 

154],    the    results   are    not   surprising.    However, 
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Figure  4-1   images  from  Figure  2-7  (top  4)  Applied  to 
Summation  Filter 
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*  ^ 


Figure  4-2   Images  from  Figure  2-11  (top  4)  Applied  to 
Summation  Filter 
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Figure  4-3   Images  from  Figure  2-14  (top  4)  Applied  to 
Summation  Filter 
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except  to  the  extent  that  blurring  is  useful  or  desirable, 
applying  the  summation  filter  to  image  signals  that  are  not 
based  on  the  application  of  the  corresponding  difference 
filter  to  that  signal  seems  to  be  of  little  utility. 

In  the  remainder  of  this  chapter  we  consider  applica- 
tion of  the  summation  filter  to  regenerate  actual  image 
textures . 

B.   AUTOREGRESSIVE  FILTER  PARAMETER  ESTIMATION  PROCEDURES 

The  first  step  in  testing  the  ARIMA  model  is  to  estimate 
the  autoregressive,  quarter  plane  filter  parameters  required 
to  model  the  real  image  textures  and  the  signal  resulting 
from  application  of  the  Laplacian  operator  to  those  images. 
For  a  zero-mean  signal,  these  model  parameters  are  found  by 
solving  a  set  of  Normal  equations.  In  these  equations  the 
white  noise  covariance  is  referred  to  as  the  prediction 
error  covariance.   The  Normal  equations  can  be  expressed  as 


^0 

s 

^1 

• 
• 

= 

0 

• 
• 

• 

^p-1 

• 

0 

(4.1) 


where  the  R  matrix  is  the  correlation  matrix  for  the  signal 
(block  Toeplitz  with  Toeplitz  blocks) ,  the  a  vector  consists 
of  appropriately  ordered  filter  coefficient  vectors,  and  s 
is  a  vector  containing  the  prediction   error   covariance   as 
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the  first  and  only  nonzero  element.    Here  a^  =  [aj^g   ^ii 
^i2  '    •    '    ^i   Q-ll'^  and  ^  =  [a  ^  0  0  .  .  .  0]"^. 

Calculating  the  correlation  matrix  and  prediction  error 
covariance  from  the  image  signals,  and  solving  (4.1)  for  the 
"a~ vectors,  provides  all  the  parameters  needed  for  the  2-D  AR 
model.  The  multichannel  form  of  the  Levinson  recursion  can 
be  used  to  solve  these  equations  more  efficiently  [Ref.  2:p. 
454]  . 

C.   APPLICATION  TO  REAL  IMAGE  TEXTURES 

Actual  image  textures  used  here  are  from  the  image  data 
base  at  the  University  of  Southern  California's  Signal  and 
Image  Processing  Institute  [Ref.  9:pp.  13-14].  The  images 
selected  from  this  data  base  are  contained  in  a  book  by 
Brodatz  [Ref.  10].  Portions  of  the  images  of  size  128  x  128 
pixels  were  obtained  and  used  as  a  basis  for  processing. 
Filter  coefficients  were  calculated  for  the  real  image 
textures  shown  in  Figures  4-4  and  4-5.  Various  filter  sizes 
were  tried  to  determine  which  yielded  the  best  results  in 
generating  a  particular  image,  and  a  quarter-plane  filter 
size  of  4x4  was  selected.  Results  of  autoregressive 
filtering  of  white  noise  using  the  appropriate  calculated 
coefficients  to  model  each  texture  are  given  in  Figures  4-6 
and  4-7.  Images  generated  by  applying  the  Laplacian  differ- 
ence operator  to  the  real  images  are  shown  in  Figures  4-8 
and  4-9.   Autoregressive  generation  of  these  images  using 
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Figure  4-4   Actual  Image  Textures 
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Figure  4-5  Actual  Image  Textures 
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Figure  4 -6   Image  Textures  Generated  Using  an  AR  Model 
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Figure  4 -7   Image  Textures  Generated  Using  an  AR  Model 
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Figure  4-8   Actual  Images  After  Laplacian  Filtering 
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Sand 


Gravel 
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Figure  4-9   Actual  Images  After  Laplacian  Filtering 
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filters  with  the  corresponding  calculated  coefficients  are 
given  in  Figures  4-10  and  4-11.  Finally,  the  application  of 
the  signal  represented  by  the  images  in  Figures  4-10  and  4- 
11  (without  the  0-255  scaling  reflected  in  these  figures)  to 
the  21 X  21  inverse  filter  described  in  Chapter  III  yields 
the  images  shown  in  Figures  4-12  and  4-13.  Comparison  of 
all  of  the  above  results  yields  the  following  observations: 

1)  Autoregressive  modeling  of  the  water,  grass  and  sand 
textures  yielded  good  results.  Some  of  the  other  tex- 
tures with  more  structure  and  sharp  local  variations 
were  not  reproduced  well. 

2)  Autoregressive  reproduction  of  images  created  after 
application  of  the  difference  operator,  with  the 
exception  of  the  water  image,  yielded  generally  poor 
results.  As  observed  in  Chapter  III,  the  application 
of  the  difference  operator  produces  a  seemingly  more 
stationary  result,  but  the  edge  structure  that  remain- 
ed in  most  of  the  images  after  application  of  the  dif- 
ference operator  was  in  general  not  reproducible  using 
a  purely  AR  model. 

3)  Application  of  the  inverse  filter  to  the  image  signal 
generated  by  AR  model  reproduction  of  the  difference 
operator  results  yielded  smoothed  versions  of  those 
results.  This  is  similar  to  what  was  observed  in  Sec- 
tion A  of  this  chapter  when  images  were  applied  to  the 
inverse  filter  that  were  not  based  on  the  specific 
data  generated  by  the  difference  operator. 

As  a  final  test  of  the  ARIMA  model,  a  64  x64  contrast 

enhanced   aerial   photograph   of   trees,   with   smoother 

variations  and  in  general  less  edge  structure  than  the  other 

images  tested,  was  tried.   The  results  are  shown  in  Figure 

4-14.   Though  this  image  seemed  somewhat  better  adapted  to 

the  model,   overall  the  same  observations  outlined  above 

apply. 
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Figure  4-10   Laplacian  Filtered  Image  Textures 
Generated  Using  an  AR  Model 
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Figure  4-11   Laplacian  Filtered  Image  Textures 
Generated  Using  an  AR  Model 
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Figure  4-12   Image  Textures  Generated  U 
ARIMA  Model 
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Figure  4-13   Image  Textures  Generated  Using  an 
ARIMA  Model 
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Figure  4-14   Final  Test  of  ARIMA  Model  on  Contrast 
Enhanced  Trees  (Magnification  X2) 
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D .   SUMMARY 

The  effectiveness  of  AR  reproduction  of  image  data  using 
white  noise  input  and  filter  coefficients  calculated  based 
on  the  statistics  of  the  image  signal  is  highly  dependent  on 
the  nature  of  the  image  data.  The  water  image,  for  example, 
with  its  smoothly  varying  and  rather  homogeneous  nature,  was 
quite  well  adapted  to  AR  reproduction.  Other  images  with 
more  structure,  abrupt  variations,  and  more  non-homogeneous 
characteristics,  were  not  autoregressively  reproducible  to 
any  great  extent. 

Using  the  ARIMA  model ,  it  seems  that  the  operation  of 
the  inverse  filter  is  very  sensitive  to  the  nature  of  the 
input  data.  Input  data  that  are  strictly  based  on  the 
difference  operator  output  can  reproduce  the  original  image, 
as  was  found  in  Chapter  III.  However,  the  AR  model  used  to 
generate  the  inverse  filter  input  (based  on  the  statistics 
of  image  signal  produced  using  the  difference  operator)  does 
not  generate  image  data  accurately  enough  to  reproduce 
images  that  resemble  the  real  images  tested. 
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V.   CONCLUSIONS 

This  thesis  sought  to  explore  experimentally  and  to 
understand  how  linear  filtering  models  could  be  used  to 
generate  texture  in  images.  Of  particular  interest  was  the 
investigation  of  2-D  ARIMA  models  to  see  if  they  might  be  of 
any  utility  in  this  effort.  Some  time  was  spent  exploring 
separable  2-D  models  to  understand  how  transfer  function 
pole  placement  affected  image  texture  characteristics. 
Image  textures  generated  using  these  models  and  applied  to 
the  summation  filter  yielded  blurred  or  smooth  textures  with 
seemingly  little  variety  or  utility.  The  ultimate  test  of 
the  model  was  its  ability  to  reproduce  actual  image 
textures.  The  purely  AR  portion  of  the  model  reproduced  a 
few  types  of  actual  textures  well.  However,  the  full  ARIMA 
model  failed  to  generate  image  textures  that  resembled  the 
source  images  used.  Many  of  the  textures  had  strong  edge 
differences  that  were  not  accurately  reproducible  by  the  AR 
model.  Also,  the  summation  filter  developed  seemed  very 
sensitive  to  deviations  in  image  signal  data  from  that 
generated  by  application  of  the  difference  operator;  that 
is,  the  procedure  seemed  not  to  be  "robust." 

Since  many  of  the  images  tested  here  have  definite  edge 
structure,  the  difference  image  had  lines  which  were  not 
reproduced   well   by  the  AR  model.   Correspondingly  the 
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integrated  AR  model  did  not  reproduce  the  original  image. 
For  this  type  of  image,  a  combination  of  a  line  point 
process  model  [Ref.  11]  with  the  integrator,  would  possibly 
have  been  more  suitable.  The  image  of  trees  had  not  such 
edge  structure  and  produced  somewhat  better  results. 
Further  experimentation  with  images  of  this  type  and  the 
ARIMA  model  would  perhaps  be  appropriate. 
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APPENDIX  A 
COMPUTER  PROGRAMS.  SUBROUTINES.  AND  FUNCTIONS 

Listed  below  are  the  names,  associated  computer  systems, 
and  functions  of  the  various  computer  algorithms  used  to 
accomplish  this  thesis  research.  All  programs  were  written 
by  the  thesis  author  unless  otherwise  noted.  Program  source 
codes  are  given  at  the  end  of  this  appendix  (except  for  the 
MAKFIL*  series) . 

A.   PROGRAMS 

1.   AUTOREG  (VAX/ VMS  FORTRAN) 

The  program  did  the  following: 

1)  Generated  a  128  x  128  zero  mean  white  noise  matrix 
using  subroutine  PGAUSS. 

2)  Multiplied  the  white  noise  by  the  appropriate  image 
data  standard  deviation  when  necessary. 

3)  Converted  that  matrix  into  a  displayable  image  file 
using  subroutines  SCALE  and  INTBYTE,  when  necessary. 

4)  Read  filter  parameters  from  an  input  file  into  an 
array. 

5)  Implemented  the  equation: 


1^1  Q-1 

y(n,m)   =   -   j-    j^   a^jy  (n-i,m-j )  +  w(n,m)   (A.l) 

i=0  j=0 
(i,j)7^0 


using  the  white  noise  array  and  the  filter  coefficient 
array  as  inputs. 
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6)  Converted  the  array  result  from  4)  into  a  displayable 
image  file  using  subroutines  SCALE  and  INTBYTE. 

7)  Used  subroutine  SUBINTFILE  to  create  image  data  files 
from  filter  results  for  further  processing. 

8)  Used  subroutine  NONC  to  apply  the  summation  filter  to 
image  data  when  necessary. 

2.   NONCAUSAL  (VAX/VMS  FORTRAN) 

This  program  did  the  following: 

1)  Read  filter  coefficient  values  into  an  array. 

2)  Read  image  data  from  an  input  image  file,  converted  it 
to  integer  values  using  subroutine  BYTEINT,  calculated 
the  mean  from  the  data,  and  placed  the  data  into  a 
real  array. 

3)  Implemented  the  equation: 


L     L 
y(n,m)  =11     aijx(n-i,m-j )  (A. 2) 

i=-L  j=-L 

using  the  image  data  array  and  the  filter  coefficient 
array. 

4)  Called  the  subroutine  NONC  to  implement  the  equation 
in  Step  3)  a  second  time,  when  necessary. 

5)  Converted  the  result  of  Step  3)  to  a  displayable  image 
file  using  subroutine  SCALE  and  INTBYTE. 

3.   CONV  rVAX/VMX  FORTRAN) 

This  program  performs  the  same  basic  functions  as 

NONCAUSAL,   without   having   the   capability   of   calling 

subroutine  NONC.   It  was  used  for  convenience  in  convolving 

certain  filter  structures  with  certain  test  images  directly. 
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4.  MAKFIL*  (VAX/VMS  FORTRAN) 

This  family  of  programs  was  used  to  create  various 
autoregressive  and  FIR  filter  coefficient  files,  using 
source  data  manually  entered  into  the  program. 

5.  SPEC0R2  (IBM  SYSTEM/370  3033  VS  FORTRAN  1.4.1) 

This  program  implemented  the  equations  derived  in 
Appendix  B  and  created  data  files  used  in  developing  the 
corresponding  graphs. 

6.  SPEC0R3  fIBM  SYSTEM/370  3033  VS  FORTRAN  1.4.1) 

This  program  implemented  the  equations  derived  in 
Appendix  D  and  created  data  files  used  in  developing  the 
corresponding  graphs. 

7.  SPEC0R3A  (IBM  SYSTEM/370  3033  VS  FORTRAN  1.4.1) 
This  program  implemented  the  equations  derived  in 

Appendix   F   and   created   data   files   used   in   developing 
corresponding  graphs. 

8.  VARIMGS  (VAX/VMS  FORTRAN) 

This  program  was  used  to  display  image  data  files  on 
the  COMTAL  (not  written  by  author) . 

9.  PIECE  (VAX/ VMS  FORTRAN) 

This  program  was  used  to  make  12  8  x  12  8  image  data 
files  from  larger  image  data  files. 
10.   INTFILE  (VAX/ VMS  FORTRAN) 

This  program  created  appropriately  formatted  integer 
files  from  input  image  data  for  further  processing. 
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11.  TRANS  r VAX/ VMS  FORTRAN) 

This  program  changed  the  format  of  filter  coeffi- 
cient data  files  into  a  form  readable  by  the  image 
processing  programs. 

12.  NSHP  (VAX/VMS  FORTRAN) 

This  program  was  used  to  convert  quarter-plane  auto- 
regressive  filter  coefficient  data  to  non-symmetric  half- 
plane  autoregressive  filter  coefficient  data  based  on  the 
transformation  outlined  in  Chapter  II,  Section  D. 

B.   SUBROUTINES 

1.  PGAUSS  (VAX/VMS  FORTRAN) 

This  subroutine,  written  by  C.W.  Therrien,  was  used 
to  generate  zero  mean,  unit  variance  white  noise  using  RAN 
(a  random  number  generator  function)  SQRT,  COS,  and  SIN 
FORTRAN  functions. 

2.  SCALE  (VAX/ VMS  FORTRAN) 

This  subroutine  takes  an  image  data  array  and 
converts  it  to  an  integer  array  with  values  between  an  input 
maximum  (MAX)  and  minimum  (MIN)  using  the  following  scaling 
formula: 


I(i,j)   =   (A(i.j)-LOW)  x(MAX-MIN)  ^  ^^^3^ 

^  '-"  HIGH-LOW  ^    ' 


A(i,j)  is  the  input  image  data  array,  I(i,j)  is  the  output 
integer  array,  and  HIGH  and  LOW  are  the  high  and  low  values 
of  A(i,j),   respectively  (calculated  in  this  subroutine). 
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This  is  done  to  provide  appropriate  values  for  image  files 
that  will  be  displayed  on  the  COMTAL  Vision  One/20,  since 
the  gray  scale  intensity  level  of  each  pixel  is  represented 
by  an  8-bit  word.  So  values  possible  (in  base  10)  range 
from  0  (darkest) ,  to  255  (brightest) . 

3.  INTBYTE  AND  BYTEINT  (VAX/VMS  FORTRAN) 

These  subroutines  are  necessary  since  data  in  an 
image  file  are  stored  in  two's  complement  form.  The  related 
variable  type  in  FORTRAN  for  these  values  is  BYTE.  To 
process  image  data  using  FORTRAN  implementation  of  the 
appropriate  formulas,  these  byte  values  must  be  converted  to 
integer  (and  eventually  real  using  the  FLOAT  function)  form. 
Results  of  image  processing  formulas  in  real  form  must  be 
converted  to  integer  (using  the  INT  function)  and  then  byte 
form  to  be  placed  in  image  data  files.  INTBYTE  converts 
integer  type  variables  to  byte  type  variables  using  the 
following  criterion  (I  is  an  integer  and  B  is  a  byte) : 

If  I  <  127  and  I  >  0  then  B  =  I 

If  I  >  127  and  I  <  255  then  B  =  1-256 
BYTEINT   converts   byte   type   variables   to   integer   type 
variables  using  the  following  criterion: 

If  B  >  -128  and  B  <  0  then  I  =  B+256 

If  B  >  0  and  B  <  127  then  I  =  B 

4.  SUBINTFILE  (VAX/VMS  FORTRAN) 

This  subroutine  performed  the  same  function  as 
INTFILE,  but  could  be  called  by  a  program  to  operate  on 


processed  image  data  arrays,  rather  than  just  image  file 
data  inputs. 

5.   NONC  fVAX/VMX  FORTRAN) 

This  subroutine  performs  essentially  the  same 
functions  as  NONCAUSAL,  except  that  it  can  be  called  by  a 
program  to  operate  on  an  image  data  array. 

C.   APL  FUNCTIONS 

The  APL  systems  on  the  IBM  System/370  3  03  3  and  VAX/UNIX 
were  used  for  matrix  manipulations  and  operations,  for 
graphing  filter  structures  and  convolution  results,  and  for 
calculating  autoregressive  filter  coefficients  from  image 
data.  All  APL  functions  except  MAKMAT  were  written  by  C.W. 
Therrien. 

1.  MAKMAT  (IBM) 

This  function  was  used  to  create  the  large  coeffi- 
cient matrices  (A)  used  in  calculating  the  FIR  filter 
coefficients  as  outlined  in  Chapter  III. 

2.  CC2  (IBM) 

This  function  was  used  to  circularly  convolve  the 
Laplacian  FIR  filter  and  its  various  inverses.  Appropriate 
zero-padding  of  these  filters  makes  the  resulting  circular 
convolution  equivalent  to  linear  convolution  [Ref.  2: pp.  70- 
72],  which  was  the  desired  operation. 

3.  GETDATA  (VAX/UNIX) 

This  function  is  used  to  transfer  image  data  files 
from  a  UNIX  subdirectory  to  an  APL  workspace. 
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4.  PUTDATA  (VAX/UNIX) 

This  function  is  used  to  transfer  filter  coefficient 
data  files  from  an  APL  workspace  to  a  UNIX  subdirectory. 

5.  MEAN  (VAX/UNIX) 

This  function  is  used  to  calculate  the  mean  of  an 
image  data  file  for  use  in  the  APL  function  COVF. 

6.  COVF  fVAX/UNIX) 

This  function  is  used  to  calculate  terms  in  the  2-D 
covariance  function  for  use  in  the  APL  function  CORR. 

7.  CORR  fVAX/UNIX) 

This  function  is  used  to  estimate  the  2-D  covariance 
function  of  the  image  data. 

8.  MVLEV  fVAX/UNIX) 

This  function  is  used  in  APL  function  FF2DLEV  to 
calculate  necessary  parameters  for  the  2-D  Levinson 
recursion  from  the  covariance  function  of  the  image  data. 

9.  FF2DLEV  (VAX/UNIX) 

This  function  performs  the  2-D  Levinson  recursion  to 
solve  for  the  filter  coefficient  vector. 
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AUTOREG 


c     This  pposra-*  ge'j?«4te:s  a'i   i^^af.e   te<th^  'jsi"j~  -^HrrE  "Oise  *s  a-j  input 

c  TO    A\j    ftUT03Er,«ESSI  VE    FILTER    .>HnSE     oasflo^ETE^S     »RE    OfaTfllNFD    F^O^*    THE 

c    filE  fil:oef,  s'J^rtjtime  pjahss  is  jset  to  generate  the  IMPJT  «mite 
c    vjorsE  IN")  su'^»:)ur ines  sche  ano  intryte  «re  used  n  prepare  image 

c     :)aTA  A?RaTS  F05  TISPLAv.  SjBOOuriNE  SU^InTcIlE  IS  JSEO  IF  AN 

c    integE"?  file  resjlt  rs  desired,  ano  su'^9outine  nonc  is  used  if  a 

c     lA'LACIA'^  inverse  FILTERIN3  STEP  15  NEEDED. 

c 

c  define  variables 

byte  aO:  127), bim(0: 127,0:  127) 

ioteqep  n,see:l»PSize»csiis«o'"l,aml,i,j,iDt,ro«,col»<ntea(0:l27,0:i 
•27) 

r»al*'<  waM,vaI2»-"(0:  127, 0:127), jsu-i<,»^(0:i27,0:l27),coef(0:o, 0:9 
!),var,msK,flin,i-nout(0:l?7, 0:127) 
C   OPEN  FILiS 

03e'>(jnit  =  l,rismp='(rarnna-tn,dataJ*'^.3at',tv3e='nei<',acce5s='direct 
;  '  ,reC3r-jSi2'»=52,'n;»Kre;  =  l2'^) 

OD*■>(J^it:2,^^w•=•(pat"■"a-^|^•^■«Mi-^al*'a.3ar',tvop=•ne«<',access  =  'di 
lPeCt',r»cor-is)2e=52,ma«r»:  =  l2''> 

03ei(5,file='  (path'nano.TataMcdbu3.:lar',status='o1d') 

C         DEFI^JE     P4RaM£TE=iS 

s«e:J=I2?aS67 
rsi ?e=l ?7 
csi re:l27 
OT«l  =3 

maM:2?<>.0d0 
•n  1  n  =  ?7  .  Otin 
c       CREATE     «HITE    NOISE     ARRlY 
Od     10     i=0,ps'ze 

do     ?0     j=0,csi??-l 
i3l:i ♦! 

caM     oqi^u  »  s  (  se;-!  ,  va  1  1  ,  V  a  I  2  ) 
Ki(«»i)=va'l 
«'>(i»jol)=val2 
20  C3ntt->ue 

1   0       C  DO t    I   TU« 

C  SCALE     «R?Ay     AN    A\0    C^JvE^T     TD    BYTE    FDR« 

c  call     srale(«.-?,»nt*3,'"'^««'n'n) 

C  call     t  1 1  oy  t  e  (  1 ->t  ?•?  ,  D  1  T!  ■> 

c  «RITE     THE     •(HITE     NOISE     I'^AGE     aRRay     TO    A    FILE 

c  do    30    i-o , rs\ ze 

c  JO    '*!)     i=0,rsize 

C  i  (  i  )  =0»  'nC  i  ,  j  ) 

C  UO          C3">tinae 

C  write(l*i«l)      ( A (n  )  , nrO ,  PS  1  ze ) 

c  3o    ro-it'p>j<r 

c  ^EaD    FiLTE-    PAJi'ETE^S    liTD    ^'■1    io-^aY 
d  3     50     i  :  0  ,  o  •"  1 
d  5      t>  0      j  =  0  ,  T  T.  1 

p?a-<(  3,  5S1     ro-f  (  1  ,  j  ) 
55  topiiar  (  T2n.  1  2) 

oO  COPit'-njr? 

50    c  D^  r  1  ->ur 
c       v'j._rTPLY    A-iITE    *OISE     3Y    ^EiL     !''A^E     3l4-;0aRT    DEvIaTID'j 
yaP  =  S3Pf  (caeM  0,  0)  ) 
30     70     t :0 , PS  »  ze 
do    90     j  =0, cs  i  ze 

•<-<{j,j)=j-«(),i'*vaP 
flO  c  ont  » lue 
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AUTOREG  ( CONT . ) 


70    c  Df>t  '  ""u? 

iPoLY     /khM?    ?J0IS£     to    the     A  J  T  3Rf  GPE  SS  I  V  E    FILTES 
ao    n  ^    -(  =  1^ ,  rs  1  ?e 
o  3    1  2  n    ••  =  0  r  r  s  i  z  e 
osuK^O  .  0<30 
3?    1  5n    i  =n,  r>m  1 

O^      1  U  (I       j  =  I'  ,  nT  1 

i  ♦  (  (  1  .ea.O  )  .ar>a.  (  j  .eQ.O)  )     33    to     lo" 

r  :>»  =  n-  1 

Col =*- i 

i ♦( (ro-, 1 t .0) .or. (col . 1 t .0) 1     30     to     1^0 

ast»'»'=(C3efCi»j)*j^(ro*»col))*risu'n 
1 aO         C  ?ot I noe 
130      cant  t  nue 

i  ii(">»ii)  =  (-1  .OdO^nsuT)  ♦«<n(n>T) 

1  ?0    Com  i  i^ue 
110  Cant  i  nuf 
FIlTE"  the  I-AGE  D4T«  ikPAV  JSING  THE  LA°LACT»N  JNvERSE  F1lTE« 

Call  ->otC  (  i  •"  »  «  TDUt  ) 
5C4LE  Th£  RESULTING  1^4GE  aR54r  tNT^    CONVE=?T  TO  BrTE  FO»m 
Ca'l   SC  3  1  e  (  »  "nsut  »  I  It  e  3»  •'^a  •  »  ml  n  ) 
Call  >  nj  Dvt  »  (  1  "t  e3»  rn  T  ) 

krjte  the  :;ese»ated  i^'.age  into  a  file 

03  153  i=Ofrsi?e 

03  1 bO  i  =0 , cs 1 ze 
atil^oiwtwj) 
1  bO    c  ont  i  ->ue 

-pite(2'i*l)  ( a (n) , o=0» PS1 re) 

150  Cont  i  -MJ; 

CLDSE  'IlES 

close(u->it  =  l) 
close(uiit=?) 
c lose( 3) 

S  t  03 

end 


NONCAUSAL 


c       THIS    o^J".^!"*    GENfWATFS    A\     MiC-;     TE^Tj5c    nSji^s    «    NOnCVUSAl    FI^    FIi_TE=? 

c  •»'«iS:     5a;i-'£TE?3    AQE    UHliI-^ET    from    a    osTa    c'IlE.     The    FIlTE»    IS 

C  ;°=LIE:     Tu    a-,     image.     SJH->iMTINE5    scale,     INI^rlE.     A';0    3  r  T  F  r,  T     ARE 

C  JSE?     TD     3C?rPipr     [Mft^E     -l4rA     A9W1YS    ^'O'    OI'^PLAf.     SjdWQijTINE     UO'iC 

c  IS    U3E5     IF    k4     INFE'-^EDIATE    FJLTE'^ING    STE^     IS    DESIRED. 

c 

c       DE^I'JE    VARIABLES 

Dyts    aO:  1^7)  ,hiTif  0:  1?7,0:  137) 

i-ite^cp    n,pswe.cswe,i,j,po««rcol»i->te3(0:127, 0:127), fsizeftfs'ze/ 
•  indf«l»in-Je«2 

r.al*''    nSLi-«»iii(0:i27,0:i27),coef(-IO:iO,-lO:lO),t«su'",-nean,dji»i3(0:i 
:27, 0:127). hi3r»,;o-.trv(-l:126,-l:l?B),^,»,i.ir. 
c       OPEN    FILES 

03en(jnjf  =  l,nsm»z'  (rat*Tna'in.i-n)oarn.'d3t'»tV3e='old',access='dir 
;ect'»perorns»?e  =  J2,-na«p»c:l?H) 
C  ooentuTits^.-ia-Kps'  fp3tTn3nn,i'n)api'"6*i.jat  '»tvoe='-'>e«'»«ccess:'3ipe 

C  ;:t',peropasi?e  =  3?»-"a«p?c  =  I?8) 

03e-\(Jr<ile='  (pwtMiam.-jatalncfi  Ic33a,dat'.«tatuS='o'i') 
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NONCAUSAL    (CONT.) 


1 nde« 1 = - 1 

cs' 7e:l ?7 

^fsiz;=-l»fsi7e 
■n,x  =  255,  OaO 

■Hi  n  =  0.  OtjO 
c        REaO    FILTER    paPA^eiEPS     INTD     4M     ARRAr 
d:>    10    i  =iif  s  J  z»»  f  s  t  ze 
do    ?0     i=-nfsiz;»*size 
reai( J»25)     coe« ( i , j  ) 
25  format  (.120.12) 

20  c  on  t  i  -lue 

10    corif  i  ou? 
c       5E4D    I^A3E    n    BE    FILTEOEO    IMTU    4N    AR^AY    AND    CONVERT    TO    INTEGE=' 
da    JO    i=0,ps>z« 

r^aH  r  i  ♦!  )     (sCn)  ,n:0,  1  ?7) 
do    'JO    j:0,cs»?e 
&!■"(!  ,  i  )=^(j) 
1 0         c  on  t  t  -lue 
50    c  0" t  i iu» 

Call  3  V  t  e  i  "  t  (  3  i  T  ,  i  r»  t  e -3  ) 
C   CONVERT  TH;  INTESEP  ARRAr  INTO  A  REAL  ARRAY  AND  COMPUTE  THE  MEAN 
'nsuT'-O  •  f^rtO 
hi qn  =  -l  TOOO.O jO 
1  3-=I  0000. 0-tO 
Do  SO  i=0,rsi;e 
do  5')  j=0,csi/? 

iTi(,,j)=*l05tfiir?-3fi,i)) 
'<SUT-=i-n(i,jl»-"SL)"^ 

'*(i'^'t,i).::r.ii-1->)  r\\^--}-'(^,|) 
i  '  (  1  t(  I  ,  j  1  .  I  r  .  1  ^,  I   1  D  .=  1  -  1  1  ,  i  ) 
bO     c  of^  r  1  ->ij» 
50  conti"nj» 

■n  -  a -» =  "  S  jm  /  (  1  2  =  -  0  «  I  ?  H  .  0  ) 
-ri  ts(  *,55)'ne3n,hi  Oh,  lo< 

S5     ♦o'"''ar('     orioinal     i-naoe    djCi-    ■n«a->=  '  . -)  1  ?  .  5/  *     h  i  qn=  '  ,  d  1  2  .  5  ,  '     lo"='d 
•  12.5) 

c     sj^r'Acr    J -lE   «FAi  cR-i-i   r-ic   o^ri 

00     70     i  =0,  127 
To     -"I     j  :  J  .  1  ?  7 

iT-(,,j)  =  i'n(i,jl-''i»a-> 
^0  C  0'>f  i  lue 

7  0    c  0"  t  i  "lu  ; 
c       FI.TER     THE     I^AGE    a5Ra«' 

Oo     ?0     nziodexl»ini»«2 
do     1  0  •)     •«:  1  n-(e  «  1  ,  ' -1'^?  >  ? 
dsu'.  =  J.0d0 

do    110    i=iifsw».fs</e 
do    120     i=T''>>7e.*siz'? 
r  o>«  =  T  > 
col =K- j 

if((po«.lr.T).or.(col.lt.O))     00    to    12(1 
if((ro-.-3t.t27).or.(coI.nr.l27))     oo    to     120 
dsu-"=f'T<(ro-»col)*coeMi,i)i*:lsu'n 

t  ?0  c  0"  f  '  "^u? 

110  C  onf  i  nue 

t  rv(n, •n)  =  5Su'" 

100  c  on  t  >  oue 

OO    c  onf  '  ""u  * 
c       CALL    Sjfl-JOjTlNE     rO     I-lv£t<SE    FILTEQ    FILTERED    I^AGE    ARRAY 
C  C^'l      noncftrvji"") 
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NONCAUSAL    (CONT.) 

4D:>     THE     «1E4N    OF     THE     I'JPUr     I^ftGE     T  D    HE    FILTERED    RESULT 
do    1^5    »=0,psi?» 
da    155    j=0.csiz» 

i  ■"  (  t  »  j  )  =  t  r  y  (  i  »  i  )  ♦  n  "^  *  r> 
155         conr  i  -»ue 
I  ^5    c  ant  i  Tu» 
CHEATt     ANJ     fMAGE     DATA     FILE 

C^ll      SUDintfile(i-n) 
SCALE    THz     ?ESJlTIN3     MaGE     AR^AV    ANO    CO^VE'T     TO    BYTE    FORM 
Call     sc  s  '  e  (  '  ■"»  i  It  eq. '"a  X  »  mi  n  ) 
Call     i  n  t  by  t  e  (  i  n  t  e3»  0  >  ■!> ) 
«RTTE     THE     ^E'JESATEO    IMAGE     INTO    A    FILE 
Oo    150    i=0,ps»ze 
da    16  0    i=o,csize 
a(  i  )=Di'T>(i,  i) 
1 60         cant  i  sue 

wpite(?'i+l)     (a(n),r»  =  0,psize) 
150    c  anr  i  ->u? 
CLOSE    FKES 

close(u->it  =  1  ) 

close(u->it=2) 

c  los^rf  J) 

sr  03 

end 


CONV 


c  THjS    P^03>'i«    GENERATES    A' 

c  y^DSE    =A?A>'ETE?5    A^E    0^ 

c  A°  =  LrE5     T-0    AM    I"tGE.     S^-,„^...._    -.^^-.     . 

c  j3E0     TO    3rE  =  A3E     I^AGE    OATA    Affk^ifs    -3?    OISP'.Ar. 


U     r^sGE     TExTjPE    'JSIN:;     ^     \0-jCAUSAl    FI5    FIlTE-* 
)?TAI\iE':>    c9?-<    A    DAlt    -IlE.    the    ciltfw    is 
;j5V0jrrfJE5    S:iLt.     I\TbrTF,     A>jU    ByTEInT     A3E 
-i  4  r  A     iBtJi^c     r -|  3    nT<;Pi   ay 


OE-IvE  \/.*.PTA  =  LE3 

oyfe  i{  :: 1?7) ,oi'CO: 1?7,0: 127) 

>iteq?r  n,rsi?e»csize,i,j.pp-.cnl,int»:3(0:l?7,3:l?7),fsize»-nfs'ze 

p^al*"^  isu'n»i-«(0:i37/0:ii?7),coef(-li:il/-ll:Il),"nsu-».rT<ea".3o'^3(0:i 
:  2  7  ,  c  :  1  2  7  ) 
C^E''  FILES 

ODff^(Jn^t  =  l,^:^•"«:•(pa^^^^^a-l'^.'•n)a'■1Tl53.3^t'»rvc>e='^l-i',acC(•ss='^)^^a 
;ct'fi"fC3PTSi7e=52.i'a»r»t  =  l?'") 

0De'>(jnit=2,nam»:'(par*i'>a-<n.iT)ari'ns'?.33t'»tyoe='ne»',access='aiP« 
:ct  '  ,rrcar3Sizr=S2»-Pa«pec  =  l2'') 

03e->(i»^i  le='  (pathi^an-i.aatalcoov^l  ,-iat  '»starus='ola'  J 
DE-I'jE  pawameteps 

rs>  7-=l ?7 

cs' 7e=l ?7 

♦  si  ze  =  l  1 

PEaO    filter    papa-etEPS    PlTD    an    aoray 

Od     10      I  :m  ♦  s i  7 p ,  f s i  Z f 

da  ?0  j=-"'s»7e.  <si?? 
r-aaC  3/25  I  coef ( i , ; ) 
?5      tap,.at  (-»2'>.  12) 
20    c  Dot 1 Tue 
JO  c  ant  1 Tu? 

P£aO  I'''A".E  rn  R£  F]LTE=EO  IMJ  A -j  iP^Ar  fi\0  CO'JVE^T  TD  pJTrGER 
da  50  i:0,psi/e 

rea-)(  I'j*n   {fl(->),n:0,l?7) 
ri  a  u  0  J  r  0  ,  c  s  1  z  ; 
t)iTi(i,j)=a(i) 
«0    c  ant  1 Tue 
50  cantiTUf 

Call  ayteinrfaiTi.inteq) 
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CONV    CCONT.) 


CONJVtPT    IH:     I-^TESE^    i^SAr    INirn    a    sfal    *p9»y    avd    CO>'PuTE    the    "tEAN 
Od     50     irO.rswF 

■''Su^^=l•n(  i  ,  J  jt-nSu"* 

bO         cant)  -loe 

so     c  onr  1  ->u  ; 

T^rt-.:flSj'r/fioat((p5>;e*l  )«(csize*l  )) 

-rite<»,55)-"ean 
55     for  ,,at  (  1^0.  I  ^) 
SU^T7acr    T-iE    HF4N    rw.-M    ThE     f  <4r.E    apRat 

3d     70      i:0,rsi;e 
ciD     "1     j  =0,cs»  ze 

i'<(irj)  =  i«(i,i)--i»aT 

''O  C  ,Tn  t  I  TU^ 

70     c  DO  t  I  -lu; 

-iL^i^    '-*£     I-«Ar,E    A.?3aT 

aa'JOnrCrsi^ff  ' 

do    1  OT    •"  =  ').  :s  we 

nsu-  =  o  .  ■)i.iO 

33llO,=„fS)7e,*si7e 

P3i<rn-  1 
col zm-l 

i  f  (  (PD-.  1  t  .0)  .OP.  (col  .  1  t  .0)  )     -3:3    to     120 
i f ( (PD-.3t . 1 27 ) .OP. (col .qt . 1 27 ) )      30     to      120 
asui<:(i-n(po-»col  )«Coef(  1  >  i  )  )*-)Su'" 
120  c  30 1 1  nue 

110  conr  » i-(ue 

djTo(n,in)rasijii 
lOOcont'-iu* 
30     conri-io» 

SCaLi     TiE     3ESJLTIN3    I-iaGE:    sR?Ar    i\1    CONVE^'f     T3    3yTE    F0»*^ 

C3II      SC  3  1  »  (  lu-no  ,  i  n  t  eCJ  ) 

Call     1  n  t  0  V  t  <?  (  1  n  t  ea»  D  '  ■"  ) 
rt^lM     rn-:    3E>;£54r£D    I-I-VG;     ii^ro    a    F ILZ 
a  0     1  5  "!     i  =  0  ,  p  s  i  ;  • 
3o     lt)3     i^T.csize 
a  t  J  )  =  =  1  -n  {  1  .  j  ) 
1 bO         c  on  r I lue 

-rite(2'i»l)      ( a ( n  )  , n  =  0 ,  PS  1 z?  ) 
150    conri-iur 
CL:)SE     -1; 


C  c 


Close(u-iir  =  l) 
Close(u-iM=2) 
c 1 ose( i) 

St  03 
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SPEC0R2 


C   THIS  PROGRAM  JCLVES  EOUATIONS  FOR  SPECTRAL  CONTENT  AND  CORRELATION 
C     IH  ONE  DIRECTION  OF  A  GIVEN  AUTOREGRESSIVE  IMAGE  MODEL.   IT  WRITES 
C     THESE  RESULTS  TO  DEVICES  S  AND  6  RESPECTIVELY. 
C   DEFINE  VARIABLES 
INTEGER  I.K 

REAL'S  PI.ALPHA.THETA.M,ALPHAS.SXM.RXK(-<9:«9).Z.ALPHAK 
C   DEFINE  PARAMETERS 
PI«J.UIS92i5< 
ALPHA«0.>»1»«» 
ALPHAS«ALPHAi>»2 
C   START  X  AXIS  LOOP  AND  DEFINE  X  VALUES 
DO  10  1*0. «9 

W.(-l.0»Pn*(2.0"PI"CFLOAT(n/9».0)) 
C   DO  SPECTRAL  ANALYSIS 

SXW» 1 .0/(1.0- (2 . 0»ALPHA«COS  C  W ) ) 'ALPHAS ) 
WRITECS.1S)M.SXM 
15    FORMAT(FI0.5.1X.F10.5) 
10  CONTINUE 
C   DO  CORRELATION  FUNCTION  SOLUTION 
DO  20  K.0.<i» 

ALPHAK«ALPHA»«K 
RX<  t  K ) • ( ALPHAK/ ( 1 . 0- ALPHAS ) ) 
I«-I»K 

RXK(I)»RXK(K) 
23  CONTINUE 

DC  10  K— 4».6» 
Z»FLOATtK) 
WRITE(<.:5)r.axK(K J 
25    FORMAT(F10.0.1X.F:0.15) 
JO  CONTINUE 
STOP 
END 
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SPEC0R3 


C   THIS  PROGRAM  SOLVES  COOATIONS  FOR  SPECTRAL  CONTENT  AND  CORRELATION 
C     IN  ONE  DIRECTION  OF  A  GIVEN  AUTOREGRESSIVE  IMAGE  MODEL.   IT  WRITES 
C     THESE  RESULTS  TO  DEVICES  J  AND  «  RESPECTIVELY. 
C   DEFINE  VARIABLES 
INTEGER  l.J.M 

REAL«8  PI.  ALPHA.  THET  A.  ALPHA  J.  ALSP  A.  ALPHAS.  ALPHAS,  TW0TN.C0S2T.H.  A.  B 
■.COS2W.SXWl,SXH,K.ALPHAK.KTHETA,ALPH2K.IC?2TMP,KTMP,INTERM.RXK(-<»: 

C   DEFINE  PARAMETERS 
PI«S.i«lS92«5< 

ALPHA»0.> 

THETA.O.O«(PI/12.0) 

IFCTHETA.EQ.0.0)  THETA.O. 000001 

ALPMAJ»ALPHA>i«5 

ALJPA* ALPHA J 'ALPHA 

ALPMAS>ALPHA»«2 

ALPHA<i«ALPHA»»<: 

COS2T.COS(2.0»THETA) 

TWOTH.2.0«THETA 
C   START  X  AXIS  LOOP  AND  DEFINE  X  VALUES 
C      DO  10  I«0.f9 

C        W.(-1.0»PI)«(2.0"PI»(FLOATtn/9».0)) 
C   DO  SPECTRAL  ANALYSIS 
C         A.THETA-H 
C         B'THETA'W 
C        COS2W.COSC2.0"W) 

C         SXWl "ALPHAS- CCOSf A )"ALJPA)»(ALPHAS»C0S2T)-(C03(B)"AL5PA)*( ALPHAS 
C     ■"CC32W) 
C         SXW»l.O/(  1.0*ALPHA<i»(2.0"SXWl  )) 

c  ws:tec5.15)w.sxw 

C   15    FORMATIFIO.S.IX.FIO.S) 

C   10  CONTINUE  ■ 

C   DO  CORRELATION  FUNCTION  SOLUTION 

DO  :o  j«o.<j? 

X.FLOAT(J) 
ALPMAK«ALPMA»«J 
KTHETA.K'THETA 
KPrTMP.t (2.0»K)"THETA)-PI 
KTMP.KTHETA-PI 

INTERM.(  COS  (KP2TMP)-(ALPHAS«C0S(KTMP)))/(  1.0*  ALPHAS- (  2.  0"»ALPHA 
■S"C0S(TW0TH) )) 

RXKtJ).(ALPHAK/(2.0«(SIH(TH£TA)»«2J))»((COS(KTHETA)/(1.0-ALPHA 
"SI )«INTERM) 
CO  TO  S6 
C   IS      RXK(J)«(ALPHAK«(2.0-Kn/t2.0»(  (1.0-AI.'HAS)»"2n 
li  M.-1»J 

RXK(M)-ax<(J J 
3  0  CONTINUE 

00  «0  J>-49,<19 
Z=FLOAT(J! 
w;R:Tc(«.<iS)r.RXKCJ) 
65    FORMaTtFlC. I . 1X.FS0.15) 
<i0  CONTINUE 
STOP 
END 
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SPEC0R3A 


C   THIS  PHOCRAM  SCXVES  EQUATIONS  FOR  SPECTRAL  CONTEMT  AND  COBREUATION 
C     IN  ONE  DIRECTION  Of  A  GIVEN  AUTOREGRESSIVE  IHACE  MODEL.   IT  MRITES 
C     THESE  RESU.TS  TO  DEVICES  J  AND  «  RESPECTIVELY. 
C   DEFINE  VARIABLES 
INTEGER  I.J.M 

REAL'S  PI.A.B.AB.AH8S.AS.BS.AK1.AICBK1.BK2>INTERH.RXK(-<9:<»).Z.A 
■K.KPZ.KPI.OPAS.OMAS 
C   DEFINE  PARAMETERS 
PI«S.l«1592iS< 
A»-0.» 
B«-0.8»»»> 
AB>A«B 

AMBS»(A-8)""2 
AS«A»"2 
BS«B»"2 
C   START  X  AXIS  LOOP  AND  DEFINE  X  VALUES 
00  10  1*0.9^ 

W.(-1.0»PI)«(2.0»PI«(FLOAT(n/»>.0)) 
C   00  SPECTRAL  ANALYSIS 
COS;W«COS{2.0"W) 

SXHl.l.0-(2.0"(A*B»(A3«B)»(A«BS))''COS(U))»(2.0»A'iB«COS2H) 
SXW.1.0/(SX>ll»AS*t2.0«A«B)«3S^(AS«BS)) 
HRITE(3.15)W.SXW 
15    FORMAT(F10.5.1X,F10.5) 
10  CONTINL'E 
C   SO  CORRr>.ATION  FUNCTION  SOLUTION 

DO  :o  j«o.<i» 

IFCA.E3.3)  CO  TO  IS 

Af;i«A»»(j»i ) 
A>cr"A»»(  j»2) 

BKl«B»«(J«l) 
BK2iB""( J*2) 

IMTERM»(AK2/C1.0-AS)-(((AK1»B)«A«BK: ) !/(I.0-AB))*(BKr/(1.0-B3) 
«) 

RXK ( J ) ■ (1 . 0/ANBS ) ■ INTERM 
GO  TO  J* 
35      AK.A»»J 

KP2>FLOAT(J»2) 
KPI«FL0AT(J*1) 
OPAS«1.0«AS 
OMAS»1.0-AS 

RXK(J)»(CK?:«KP1»AK3/(2.0"OMAS5)-((FLOAT(J)«KPI«AX)/;2.0"OPAS) 
•) 
5(      M.-J-J 

RXK(M)<RXK(J) 
SO  CONTINUE 

DO  60  J«-«9.6» 
2»FL0AT(J) 
WRITE(6.6532.RXK(J) 
45    FORMAT(F10.1,1X.FJ0.15) 
6  0  CONTINUE 
STOP 
END 
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PIECE 


c       T-ijS    p^03al^'•l    T4K=;-,    4     i?fl<i3A    =T£:E    0  •!    oe    t    SlP'Sl?    image    O^ia    FILC 
c  a\iO    r'_A:£s    THiT    Oili     I'ltn    4-1    i-'.'GE    -IlE 

c 
c 
c 

c       DE'InE    y/ft-^fA^LFS 

o^te    sO:Sin,oim(0:Sll,0:Sll),ao{0:I?7) 
c       OPE'J    PILES 

ooe-<(jnit:2.r>4mi.=  '(pat'>Tia-iri.i'»)trees64,-lat'»tvOf='r>ew*»access='cJ»r 
lect'/rerorasireslfj/iawr-c^oo) 
c       REaJ    IV1A3E    TO    SE    FILTERED    P'TO    A-,    aRSAr 
do     50     i=0,l27 

p-a3( 1  '  i  ♦>)     (a(o) ,n  =  0, 1 ?7) 
aaJOjrO.l?? 

0  1    Tl  f    1    ,    j    )  =  .^  (    i    ) 

a 0         c  DH  r  1 -1UI" 

c       .^''IfE     THE     DESISED    I'^iGE    °0?irO\     INTO    4    FILE 
dD     15")     i=b3,l?o 
03     16-)     i  =^S,  t  ?6 

aD(«)=Di'"(i,j) 
160  cant  i  sue 

wr  I  t  e(  ^  '  1 -o?)     (  a3(n  )  , -1=  ■>,  b5  ) 

ISO     c  ^T  i  -<u? 
c       CLOSE    -I. £5 

close(uTir=l) 

sr  03 


INTFILE 

c        TMiS    P?3;RA^■•    C'-'E^TES     i     l?fi<I?S    OP    b-itb'i     I-JTE^E'    PRE    F»0-    i     129x13^ 

c  lN  =  ur     l^'iGE    ^IlE 

c 

c 

c 

c 

c       OE-I>jE    V49I4  =  LES 

Dvt;    3(n:l27),t)i'n(0:157.0:ip7) 

lot-a-r    i,i,r(,«,l,)r.teQ(0:l,?7,0::?7),suK»iiai,lo- 
r»a\»'     iiear< 
C       OPE -J    FILES 

03e->(jnit  =  l.n5-ie:'(rar'>'(a-.n.j^)-)»«3r3.iat*,tvoe:'old',acccss=*d 
;>p?rt'»recor3siz»=5?»'<a»p»c  =  l?9) 

C  3Den(2,fti»='  lrarMia-<-i..Tata)'5>"-i?.d3t',s?atus='new'> 

c       ^EaDI«i:;ET^|^EFILrESE:^INI04■.    A«::ttr 
33     50     i  =0, 1?7 

r?a3(  I  *  1  ♦ 1  )      (fl(n  )  ,n  =  C,  1  ?7) 
do     JO     j=0,l?7 

Oi'flf),  ))  =  =»()) 

il  0  C  3n  t  >  TUe 

50     C  3nt  I  -<u; 
C        CUvjVERT     THT     ayi£     *5Rnr     pj  T  3     a»V     IVTEG:^     «o^4r 
Call     3vte»nr(TiT,t'ite3) 
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INTFILE  (CONT.) 

c   ZO'^P.nZ     HISH  AND  LO^  VALJ£5  DF  Th£  I'^AP.E  DaTA  ANn  ^RI^E  TO  Tt='*'lNAL 

h  1  3  1  =  '"' 

1  o-  =  ^55 
To  51  i  rO,  la? 
rlD  5?  i  =3.  i  ?7 

i'(int?o(i»i).lt.lo.)  lox^'nteaf'.j) 

51   C  on  t  1  -iu« 

-r  1  t  «  '  *  ,  55  )h  1  3^1  »  1  .^« 

c       CO-PjT;     -EiN    3F    IMEGE^    as^flf    an^    su=!T?acT     IT    ffir\--t    THE    OaTa 

C  50"='^ 

C  -JO     55     i=0,l?7 

C  :)o    "5     )=■).  1^7 

C  su'«  =  su-n  ♦  1  n  t  ea  (  »  »  i  ) 

c  us  CO-itinue 

c  55    ;o->rinjff 

c  K.a-.=  M  5dt  (su") /•(  1  ?'^."»  12?.0) 

C  5  J      f  Or-nst  M   1  0,el  0.  5  ) 

c  33    5-5    ^  =  -:,  127 

c  no    Jb    i  =n , 1 ?7 

C  int»-i(j,j)rint(»T('«jl-int(TCsn*U.5) 

c  5r3    :  o->  t  <  -^  j» 

c  ARrli     JNilESE-     A5;}Ar     I^JTO    A     l?6<l?f»    Dili    PILE 

c  k:i2'< 

c  1  =  1  2  e 

c  ^Pi  t  ?( ?, 5^)k,  1 

C  53      for-nst(2iS) 

c  30    5:     i  =0/ 127 

C  3061      1=0,  112.  It! 

C  «r,te(2/7'"i)i-.t?-3(i,)),inreT(wi*'.  )ri-ite3(i,j*2)»<r"te::fi/i+5).i 

c  ■-,te3C',j*'J),i'^teT(i/i»5Ki'~it»a(i,it3)»'"teo(i»it7),i-te:;(i.i*6).ii 

c  •re3(>.i+^).i-ireg(i,)*lT),'ntccifi,)<-ll),)-it?q(i,i»12),irt»q(i,jtl5) 

r  •,)-!t»g(i/i*!-),>->te3(i,itlS) 

c7t  ♦of-natflaiS) 

c  60    continue 

Z  5l   ;3-iti^je 

c   CLDSE  -I-£3 

close(u-iir  =  1  ) 
close(u-iit=2) 
St  00 

e-\a 


TRANS 


THrS    P^O'lfA'^    ''Ea)5    OaTj     "IHjT     tc    a    pile:     in    P'EE    fOff'^AT     A  1,0    CONVERTS 
IT     TD     A     oROGtJA-    REaOASLE     fl-^'&J. 
i  Tt  eq-r     i  ,  j 
real  » 'i     a(  16) 

oo'*n(l#f»le='lp*th-n4n-i.3atala3-at.3it'fsr»tos='old') 
oo*-i(2r  fi  le=*  (P3t>''"»n-..-)ata)  'ciwar  .-)a'  '  ,srafus='ne-'  ) 
rea3(l»«)i,i,a(l),a(2),a(J),a(y),a(S),a(6),a(7),a(fl),a(<'),a(lO),8( 
*ll),a(t2),a(I3),a(l'i),a(15).a(l6) 
do    10    ic:l  .  Ife 

«ri t?(2.20)3(k) 
20  ^o'-'-at  (a20.!2) 

10    c  0"'  »  "<ue 

St  00 
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NSHP 
c      THIS   p?j;p!M  GeNr=sr^3   4  J    T^aGf^    rEtrior   'ISIS-    vHrT"^   NOISE   45   4 \(   Input 

C  TO     a\;     VjTO'EG^rSSI  ;£     FILTE'     vHisr     3i540-rTr:?5     jpr     OHTfirrjfr?    p:?0.«    TH^ 

c  ="Il£    FILIOE'^.    SJ^W.l  ITImE    PGi'iSG    IS    JSE:^    TT    GF'JcRarp    rnr    jrjpjT    /«Mir£ 

c  'JOtSE     iN?    SJ^^DUtlMcS    SC4LE    iN'^    IMTBrfE     »»-;      ISED    TO    PHlPA^i     I^4GE 

C  -JATA      45R4YS     F05     9IS?L4T. 

C 

c       DE-PJE    VHRM^LES 

Dvtf    i( d: 1  a?) ,dit(0: 1?7,0: 1271 

:27  )  ,  -nam! 

p?dl»8  tfall,v3ia,-i(0:l?7,0:l?n,isu'"r>'t(0:i27,0:i27),coef(-2:9,-3 
;  :  9) 
c   OPE'J  FILES 

O3e">(jriit:l,n8'"f='lpdt"''a">".c3atal»n.dat'«tvoe:'ne«',access  =  'direct 
!',reC3P1SWp=i2.m-»«rec  =  l23) 

03e->(j"it=2»n3me='{rat^iT<a-in.aat3)potat;.iat'rtvop='ne«',access=*ai 
1  fect  *  I  PrC3P-isi  ze  =  52."axpe:  =  l2P) 
03e-i(5»fite='(patn-iianr>,t(ata1fc5a.3at'.st3tijs  =  '3l:3') 
c  03?-i('J.fi1e='lP3tn'"an'>.'3atalts.Tar'.status='"e«')  ' 

C         DE'lNiE      ?4RjMrTE;3S 

s  ?  e  3  =  1  ?  J  ■J  5  6  7 
<"si  7^  =  1  ?7 
cs<  7e:l 27     - 
o^l  =2 

iijw  1  :  -  1  tOT  1 
C         C'^rE     A-«IrE     NJI3E     43S4Y 
Qo     lO     i:0,psize 

3D    ?n     )=0,csize-l 

c  "<  1  1     03aijss  f  se»  T»  V  a  1  1  ,  \/ a  I  2  ) 

«-l(l,   I   ):»/:«)   1 
-->(.,  lOl  ):v3l  2 
20         CD-'ti-(u; 

10     Cont'-iur 

c       SCiLi    s"4r     'i'i    4\i0    CTMvE'I     TD    urTE    fds^' 

Call      SCjle(-n,inc«3) 
Call     infDvte(i'it»a»DiT>) 
c       "^HE     THE     /v^-ITf     MOISE     I^^iGE    irjcar    TO    4    FjlE 
33      50     i:0,psi7e 
30     JO     i:C',cs<ze 
a(j)=oi'n(>.i) 
uO         c  D1  r  1  -lu? 

«rife{l'i»I)     (a('^).->  =  0»PSI^e) 
50    c?'>ri"nj? 
C         'ESD     FILTE?     =4?a^-£rE55     I;T3     4N     itJR.-Y 
ao     50      I  :-noT  I  /  n-n  1 
3-)     SO      f:0.3T.l 

r-a-iC  J,  ^S  1      coeM  1  ,  i  ) 
5S  f  3P"^t ( t2: .12) 

5'"'      C5it'1U» 

C        AH'aLr     'HjTr     r.  j  ;  3  E     TO     1 -lE     -  J  T  D  ^^  r  ;  ^  -  5  ,  ^  r     r  ; ;_  7  r  o 
To      11^      ••  =  ^'»CS1?* 

ao    12'^    o  =  0»psize 
oso»  =  o.Odr' 
33    IJO     i=»<o-"l,D-l 
Tj     lan     j=n,Tnl 

i  '  (  (  1  .ea.  ^  )  •  a-^d.  (  j  .^G,  .^ )  )     33    to     MO 
PO-  ^"^  - » 
CO  t  =■"-  J 

i  f  (  (  P3«.  I  t  .  0  )  .OP.  (  coi  .  I  t  .  "1  )  1     00    to     \'*0 
dsu-<={coe*  (  «  »  i  )  •!  »<(po-,col  <  )  ♦-isj-^ 
1 uO  c  0" t ' 1U? 

130  c  on  f  i  nue 

,«(t, •»)  =  (- 1.0d0«nsuTi)»-n(n,-") 

120  c  o'l  t  i  ooe 

I  1  0    c  o"t  1  ->wis 
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NSHP  (CONT.) 

SC4LE  TM^  RESULTING  IMAGE  4R?»r  ANO  CONVERT  TO  BYTE  FOO'^ 

call  sc5>-(i^,iTteq) 

Call     intcvtp(inre3/?«Tl 
.tSTir     Th:     3E'i£f?4rE'^     !"^4GE     INTO    A    FILE 

do    160     j=0,csi7e 

loO    csnriTu? 

Mritef?'i*l)     (  a  (  n  )  » -<  =  0  »  rs  1  ze  ) 
150    cant  1  iiij 
CLDSE    -IlE5 

close(u->i  t  =  1  ) 
Close(u->»t=2) 
c)ose(3) 
CloseC) 

St  03 

e-id 


PGAUSS 


SJ^ODMTINE    03iJSS(-<.Zl  ,Z2) 

I  N  T  E  3  E  ?  *  J    -^ 

A    :    ST?r(-?.non    .    4L3G(Wi.'Jf  •< )  )  ) 

3    =    t,.?851S5:)0    *    SiN(K) 

Z 1     =    a     •    C?.5{>i) 

7  7    :    i     •    S  I'J  f  B  ) 

^  £  f  J  ^  N 


SCALE 

TnTS    S  .H-'.'J  J  r  I '.E     ^C'iLtJ    f^i    i'"if     ru     r:TE3P-<    valJES    -'ET'.EE';    -     '.  I  v  E  "i 
^fl  ,  I  V.J  ^     I'J  )     '1  .'I  -"J-' 

lot^q?"     i-it"»3(}:  I?/.'"*:  1?7),  I  ,  i 

P«a1  •'     srr(':  I27f0:  '^7)  r^iiTn,  l3w»'ne/in,na«»'''in 

=  E?p:)9'    3CaLr*r. 

ni3i:-tlO.">n.OlO 
1  3-:l  OOOO.O   tO 
dD     1 0     i :n, 127 
33     ?0     j=0.1?7 

injrr(i,i).l-.lo-)      l5-=arr{i,;) 
if(5Pr((,j).n«.niqn)      •>iTn  =  ar<"fi,j) 
Sj'»i:arp(i,j)f5oTi 
?0  C  ont  >  TU? 

10     C  3nt  1  iu» 

2S     ♦or-natt'     or»srate-<     imaq»     jar-»-    ^  i  -jn  :  *  ,  j  1  2  .  5  ,  '      1  o- =  '  .  d  1  2  .  5  ,  '     n>«an  = 
I  ,  3  i  ?  .  5  ) 

CQ-jTroE    SrALllG    ftNiO    CDNvE^r     to     I^iTESE?    FD?*'' 
do     iO     ir 0,127 
dD    'JO     i=0.l27 

i-it«:j(i.i)  =  )nt('(arr(wi1-)o-**M-r<a«-''«j->)/(ni-)h-lo-)))+0.5dO) 
iTtr3(i,  i  )  =  intej(i»  i)»int('nin) 
40  c  on  t  i  -<u» 

^0    c  -»">■'  Tur 
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SCALE    (CONT.) 


50 


65 


1  1  j'>  =  -  I  lOCO  .  OtO 
I  0"=  1  JOTO.O-iO 
0-)    50     »  =0,  1  ?7 
33     'yO     ;='),l?7 

•  <fH03tdnrsT(  i.  j  )  j.^t.^noh)     S,-iKrf|,ar(inre3(i,i)) 

S  jmr  f  1  o^t  (  i  nt  •-3(  i  ,  j  )  )  tSo* 

c  3"  t  1  -Mje 

C  ant  1  ->u; 

-ri  te(  *,i>5)hi  3n,  1  D«,m«a-» 

♦-5rT,at('     scaler)     i-»3n-    Tsti-    Niq-^:' ,  .11?.S,  '     Jo 

:2.5) 

e-i-J 


».112.S»'    T<ean=',dl 


INTBYTE 


HIS    Sjd^OjMNE     TA^Es    4'.     I.r:::"    Assar     A-O    CDMvEPrs     IT     I.,TO    A 

ay  Tl  a  V«  •  Y 
:£-I\;£  VAOifi^LcS 

1-itsq-r  i,i,i-ifa(0:l?7,ri:i27),n 
ovf?  DiT,(o:i?7, 0:127) 

n:0 

a  3  10  ) r  0  ,  1  ?  7 
Oo  20  j=C, 1 ?7 

'((lit-T(i,i).lt.0).Dr.(infe-)C..i).gt.25S))  n  =  ntl 

fCmteod.iKlf.D    oi^fi.ij^-l?-* 
♦(inreT(i.i),qt.?S5)     Di-n(i,i):i27 

M(»"tea(i,iJ.le.l?7).:»o:,,(i-,teqf>»i).3*.0)J     oi'«(t,j)  = 
I  1  ->t  •3(  i  ,  j  ) 

i'((»-.t-3(i,j).3t.l?7).an-i.f(-,teqfi,j).)e.255))    bim(i,j); 
;  j  -1 1  •■!(  1  ,  »  )  -25t> 
20         c  0"  t  i  -<ue 

-rite(*,30)     n 
IQ     ft.ri.at('     THE     'JU^'SE^    0-     ^OP.'TS     DiT     Qc     i?A\G-:     IS*,i5) 
r  •  f  jr-i 
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BYTEINT 

c       THIS    o?03W4"    Takes    A    ^tT£     &n(;J^r    a  "O    :i)'"VF=?rS     IT     INTO    4N     jrJTEGE^    array 

c 

c       DE-riE    ViRMSLES 

i-^tiasr    int-:j(0:l?7,0:I?7),,,  ,,n 

0  V  t  ?    0  1  '^  (  0  :  1  2  7  ,  n  :  1  2  7 ) 

-      pEjfts*"    :0\IVE'^S  I  on 

nzO 

3o     10     1 :0, 127 
do     ?0     j=0, 1 ?7 

if(f0»-'(t,j).lt.-l23).or.(0i*(»,i).3t.l27))     n  =  n*l 
iftoiTifi.jJ.lt.-l?")     inre3(i»i):-l?^ 
»f(3i-i<(i.j).3t.l27)     i-<teq(i,j)  =  l27 

i«((Di-n(i,i).qe--l2?).anC.(niT.(i,j).)t.O))     inten(wi)  = 
:oi  •«(  1  .  j  )  ♦25o 

iM(r>»-n(),j).3S.0).3n3.(Ci-n(i,j).le.I?7))     inte5(»»i)=5«'"(i»j) 
?0         CDitiTue 

10     C  0"^  t  >  -lu? 

-  r  1  r  -  (  •  ,  J  0  )  n 

30     *Tr^at('     HE     ^JU-"?E'    0"     =OPJTS    C'JT    D?    ^AMGE     IS', IS) 
r  r  t  J  r  -> 

SUBINTFILE 


SjOrojti'^-     S 'J  D'-'t^i  !?'<•") 

c       Triis    S.'3?3iII'Ji     :S;ArE-     i     12^«;2=     ^^r■.ZE^     ^jta    rlL'^    fdt"    f'lLTE^E:) 

c 

c       DE-l-'ii    vi^ffl^LES 

i-it;a?r     1,  i,n,.,l,inte::('):127, 0:137),  5UTI 
r.,^i«-    ..e3n,)-i(0:!27»0:i27),h,3r,,i3, 
c        O^EN     rJi.ES 

03»?->Ci»file='  (rat'n<a'">.iati)'li'inar.-!ar  ',st3tus='"»-') 
c       CO-'PjTE     -iI;h    ATJ    LDa     /iLjES     D*^     T-iE     I'^AGE     OaTi     4\/0     ^^ITE     TO     TE^'^TjAl 
H  i  Qr-:-  1  .  CTb 
1:.-:  1.0 -It, 
3  0     31      1=0,127 
C3a32j:  0,127 

if(i'"(i,j).lt.l3-)      lo«=i'"M,i) 
'♦(i'n(i,i).Tt. '>>;)'>)     ^iqi=i'nM,/) 
52  cant'-iu? 

31     c  D^t  1  ->u? 

write'*, 55)bi an, Io« 
?5     <  :r.riiat  (  ?J2C.  12) 
c       :Oj='j''£     ^Esn    Qc     II<tE3E^    i'-'^isr     fiw;:    SU^IPacT     IT    rP^M    ThE    OaT.i 

c  aoi5i:T,l27 

C  io    us     j:0, 127 

C  su  n=  su-1  ♦  )  n  r  e  3  (  1  f  i  ) 

C  '^  5  C  0"»  t  i  "ue 

c  J5    ro->tinje 

c  Kra-<=flTar(su'')/(l2-.'i«l2°.  0) 

C  ,<ritf(«,3U)sjm,'ne'*n 

c  5-1     f  o  r  -n  s  t  (  i  1  0  ,  e  1  0  .  3  ) 

c  nDisi=.T,l27 

c  a  o    «J  S     i  :  •1  ,  1  2  7 

C  )nfeq(i,))  =  titea(i,J)-int(-ne-»n»0.5) 

C  "  S  co-It  »  rtxje 

c         is    ro->tinie 
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SUBINTFILE 


■Jo     'Jb     >  rO,  l?7 

d5    55     ; =0, I ?7 

'-'f?a(i,j)  =  int(i-nfi,i)»0.5d0) 
S5         c  onr  1  -lue 

^WITE     iNfLSEP    ARPAr     I'JTO    i     l?M)fi?B    yars    fliE' 

"  r  1  t  e  (  1  ,  5  =5  J  V  ,  t 
3=?     <opi.at  (?i  5) 
03     SO     1 rO, 1?7 

To    i>'"'     i=0,  112,  IS 

..  .nf  eT(  ,  .  ,  My  ),  ,ntea(  .  ,  i  »  IS)  ^t;.v,,, 

70  f  3r„at  M  oi  5) 

t>0  CsntiTue 

50     C  on  t  1  -\u; 

CLDSE    =1.E3  t 

close(u-iit  =  l) 
Pf  t  jr-i 
e-i'i 


NONC 

Sjorojtioi^     no-(C(tPv»>'«) 
C         THIS     SjSVOjIINE      T4<£S     A.g      1-43E     Ao;?-T     A-JO     flLTE^S     IT     -«ITH     A 

c  MUsCiJ5A.    fl"    -ILTE^. 

c 

c       DE'lNit    V4SIA3LE3 

real**  -)Su-n»i''(T:i?7,5:l37),cor*f-10:W#-lO:10)»'nsu-n,i>ea"»3j'^D(-l: 
:l?8,-l:l^S),ni3^,)o-,trf(^:l^7, 0:137) 
C   OPE"^  FILES 

03e-i(l«file='(ratnTiao->.:Jaralncf>lco?1.3at'»statuS='oH*) 
c   DEFINE  PASA'^ErERS 

)  1  CJ  •  «  1  =  0 

1  -i-ie«?=l  27 
rs  > ;e= I ?  ' 
c»i ;e:l?7 
♦  si ie  =  I  0 

c       9EiD    FILTER    °A,iA-ETE0S    InTD    a.j    asSay 
anlOirmfsiij.fsw; 
Od    ?0    j =" ♦ s 1 ze# f 8 » re 
rsa:J(  I  .  25)     coe'  (  <  .  i  ) 
?5  f  3P..at  (i^O.  1?) 

20         cant  t lue 
1 0    c  snt  i iu» 
c        TA<E     T-iE     I<jPjT     4?»Af     A-gO    CD''  =  iJTc     ThE     •*£**■*.     Hlin,     AND    LO"    VALUES 

h  1  gh  =  -  I  0000.  0-jO 
I o-:| 0030.0HO 
03     50     1  I  I  ri3i?  «  1  ,  i  n-(e  •  2 
03     iO      /  : I nOf  1  1  .  ) n Te « 2 
•rSui  =  t''v(  i  »  j  )  ♦•flsuw 

if(t«'v(i»i).3t«'''3^>     hi3h:trv(i,i) 
i  f  C  r  i*y  (  i  ,  i  )  .  I  t  .  I  o«  )     l3«  =  tpv(i,j) 
60  c  ant  <  ->oe 

50     C  3nt  I  TU» 

write  (*, 55)-"»»r>.'<'3"«'0- 
55  f3riiat(*  inout  i^aqe  array'-  »ean=  '  /  d  1  2  .  5  »  '  h  t  o"  =  '  »  cJ  1  2  .  5,  '  low=*,dl 
.2.5) 
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NONC  (CONT.) 


30  ''O  n  =  1  nTe  X  I  ,  1  me  «  ? 

do  100  ■»<=irne«l»inaex? 
asu«=0 . OiO 
dD  IlO  I='^♦s^^•.♦sl^e 

r  o«  :•>  -  i 

C  D  I  ii- j 

i  ♦  (  (  r  D- .  1  f  .  »  n-le«  I  )  .  or  .  (  CO  1  .  1  t  .  t  nd?«  1  )  )  OO  to  1^0 

i  ♦  (  (  row.  q»  .  I  me«^  )  .  or  .  (  CO  I  .  3t  .  »  nl- «^)  )  30  to  1^0 

dSU-=ftrv(PD«<»Col)*COeMl»i))»3SU« 

1 ?0        c  0" t  »  nue 
IlO      continue 

i  ■«  i  •>  r  -n  )  -Is^'" 
100  contiiue 

90     C  on  t  «  "»u« 
C        COMVE^r     THE     t5O)fl50     iR^At     INTO     A     l?8ifl?8     A9S4Y,      JF     NECESSARY 
C  do     IJO     »:0,l?7 

C  30     lyO     j=0,l?7 

C  1  •«  (  1  »  i  1  =  3U-KD  (  i  .  j  ) 

c       1«0         roTtiriue 

CliOcO-itinje 

c 1 ose (  1  ) 

return 

end 
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MAKMAT 


[11] 

[12: 

ri3] 


P«-nAknAt.niC>;  A, 


nil.- 

A--     *  1 
A  r  ')  ;  0 

ato.  1 

'  t  •  1 

TTARTZ'^'I 


c .  D .  e:  .  F  ,  G  ,  H  ,  I  ,  J  , 


'- .  n ,  0 


'  1     66 
1>-4 


LOOP 

3: 


;"TARTlt-5TAfiTi  *-(I-1 
SJhRT2'-STpiRT2->-(.1  +  \ 
R0U'-"1 
LOOPr  ■  ROU'-ROU+I 
-«vROUi  =  ';)  '  zROUi 
-t(  ROU=  I  )  /ROUI 
COLi'-COLi*i 

cdl;'^<:oli+<  i-«  > 

CuL-*«-CuL3  +  1 
CgLS^-rCL'^t  1 


1  > 

) 


n- 


100 


CC2 

^    J*-<    CC2    H.  I  1  .  I2.3I0.K'i  ,N2.  Ti  .  72 

n  ]     fl 

r2j       «    FUNCTIOr)    TO    DO    2-D    CTRCULi^'F:    CONVOLUTION 


f?' 


'4l  -k(   I   »/(  cX  ^.-;oH;  -2=0  c  <  I  /r';:'2  I  fj 

m  D- '  InPROPEr.  ARGUhLNTw" 

CM  -u 

r?]  PL":: IN  nio*"^ 

["1  fj'  «•  1  to  . 

F'^l  N2*-l-*'gy 

[101  rv  <f)j  ,  Hi  ,1  CO 

Mi:;  Jl'-O 

ri2l  LOOP'  Ti«-(- 1 II -M  Mij)a)H 

n^j  ii-o 

T'  4  I  LOOP''  T2'-<'-'  12*1  )  )eeTi 

r<f]  Yri?;  11  i«--/  +  /yx  r? 

fi  '_.1  -»tN2^  i:-^i2»  1  '  i.no--;' 

1171  ->■  f'l  M  1  '-I  1  *i  '  /Ltujf-' 


GETDATA 


■    --'i._    L -^         ■•■--.,--       -:.■  ..I..        f-'r.-f.  ...:n. 

»i    _  _ 

rr       f  L''~-f  T"     Hi'*':'  ^-i"  ':^^-i"'i:'-  ^~   r.-_   _"l':     r?"-:'-,"     Eriior;     i".. 
:r;     f>  rj"K'.;T  f.ecoM'  ij  i'hafe  ct    ti-,c   u--^ 

rj^rE-^'"^!'^^ ,'.'-''     '  *NAr?t  '  ■'  '     Dr.  '^A  ' 
^.  r  .  -  A.I  >i  ri  r_     '     (  3  ~  ''■ 

;  *  ■' ;     -ij  •":  :^%;:'-":  11^  ZJJVC'   'Rir' 

'"2"       "^  ■- .'"YJ"'!.""^  ""T-'      'ETF.'^.E '.  7~ 

"■'"■"'    --lC-'.,-' ^PZL  ;  !•• '-^  ■  J/.'L'.TEI 
,-,-=.       L.I"£'  ■  |_iNE- '  -     •  '  i 'JL-r'~T'- ■  ~  ■ 

'22"       -"-■"O"- 

"  2 " '  ■■    "  L  u  -T  _  r ','  r ,".  c :"  r  ■  -  r  ^'  v  r^     p  e  •" ' 

"2^^    E^I      F<-Nno^ 

Z2ii  EPkEYIT  'FILE  DOET  not  EyZr"^' 

PUTDATA 

'    r    p'JTtf^^A    NArf    ZIO  .C  .  0"i  ::kESf     fE-'.NPEC.I 
*  ~       p  _    _  _ 

2'       fl    ^L'f'2~''0"    '"0    UF'"''"    PAT^     .VF""''-  ,     A.PP-  ■       E'C.  ■     IN    -PlE    FOPi-iA       C.'j    T.iJ    DITl 
2         fi 

-  "  p  r :    ..  . 

■"  "i  fi_:-- ^(M^(-,c     .  .,  I  yo.^rr  ■  '  '     L'.-.T''  ' 

C'-  '  2fi;."  ' 

"':      rvr'-E--'- '00  2rv''i    2' 

-. -;         f ..    epaTE         ."J-^rL 
■  '         :'2;'PT  :f-'t''y  '''^'    ' ' 

"  -  .'i  'l  r.  ;   I"   -Uf  n  "'T.""  ' 

■  *  -i'T     i" -'I'" '    ' '  '^"''  _ 

'■  :     '       C'^iz^C " '^  (  >  /  -jF     •  -i  

^  .  2  ;       i;.-~''..4  I  .'  i  .  .  „F  ,  .-•. 
r  .  cr  -         r  -  •  ''PEC  ,  -^  .'  or 

:■.-!    i_rM-.r.     MfJREC  'I'  1-1  '/r"  !irH 
:if':      PE'2-;     MFlil.T 

•"  -;  I      :":'PE^P-Oi'VP    'C'E:  ' 
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f 


COVF 

v'k*-SlZt   CCVF    F.DI0,'i,l-.L,K,R1 
FM       «    t^UNCTION    TC    GENERATE    TERn^    IN    2D    COVARIANCE    FCN   FOR    IHAGE 
4]         F«-F-MEAN    F 

6]      L<-rizEroi 

7l    fC«-.7IZEhl 
S'i         P'-fi-l^-riZEpO 

1V1  L00^-2:RrLiKD«-(*/  +  /((L,K>+F>x((-(L.K)>+F>)-ic/pF 
I2I   R1  ^L,K>^♦/■^/<  (L,-K)4.F)x<:(-(L,-K)  )+F))-x/pF 
I3I   -»(K>K«-K*1  )/L00P2 

.14]  -kc>l*-l*o/loopi 

CiSl  fi*-lP    0  1  *Rt  )  .R 


CORK 


n    Fi  (s'TZCr    '■':■    F.r~lr>-'~.    il    CC'J,' ■'-    '-.<'<:Z    F'jr":~  ZC-"    -' 
_.1   » 

■-.-■'   niL— :> 

1 1  '•   Hr--;izE  'j.'v-  c 

»■!-:, J ZtLi'-'  • 
.■■  "i  '■r-siz:''  1 

"'")  L'-O 

.',,1  LOr,c.o.c.;r|-|_  _,,  T^^_K.)f  <-K)<»R[L.  1 
"121   -<K>K-K*r-/Lur'P2~ 
''13D   -•(L>L<-L-1  ?/L'JDP'< 


r  1  T   fl 


MVLEV 

^  AH^^ni'iliV  RFT,RT;N:N&;NDrh.  AA1  ,AF;  AB:  AF1  .  AP1  ,  EF  .  EB  .  GAhF  .  GAMP  ;  DELT,  I  .DIG 

n  P"^-"^  I"  THE  fhA'^PIV)  rORR'ELAT'ON  "^UNCTION.  PlAME  J  CONll^'^r    p'.  j). 

fl  ~C  'i^P'TArf'  'Cr'^-r'Rl'    PApMrETtRi"    ,--,r ---..'■'   ;1,  PACI'-UAP:.  Ahi-nAi'  ^ , ,  .z, 

n  EPROP  CGVriP;! ,'-'CZ  E'^-«'^rM.  1.    PU~  -FTA..!'-!  TO  QPTAIN  lN'"JVATiCr«'i:' 

«  ^ilTER.  lZKEUITE  FOk  PACKUAPD  f»Rf-nc:~FAS. 


1-41 

rr""  aiO'-o 

rf?T  -.  I  \<  f  opp-^>roT;^7='>oR-'T>/rTART 
r9-i  'INVALID    ARGUrt-IT' 

r  in ";  -rt 

hi « T  H 

C<2l  rTART.N'--<-(oRPT)COD 

I":  3^  Nf(-<,  oi--p-)m 

r<  u1  F:Tr-   0    2    <  "«}RPT 

.'<  5  J  I-'  >  NP;  •  .  =  1  fJp 

[16.  -<N>1  ;/KE:'JRr 

ri  71  hh  i«-'- 1  .nb,nb/pI 

r  •■;••;  E--?p-^ro;,j 

,-201  E.^^-E- 

"2*';  -rGHALG 

[2^1  REC  ..P.:    AAi«-riVLEV(N.Nt.Nr.-;tR'''T 

:"  Ji"  Zi">-Mr<  i  1* 

:"■*  -■-   r^v;:-\ 

'  ::i  AMI  ro.  :  ^^-f 

:  ■••^.■  J  r-r,  ; ",,  .     -JO  i  fi"-Hxf)£i 

["-•'J'l  D'^'  Tr  <  i^i  Hr..:r. .  NF- 1  "  ■  1    0   f'    -tRPT >.•-.' v'IDIn.  NF>;  poAF* 

r"--  1  r.r-rr-.  ^FI  '-.xOELT 

r32l  GAMl<-^EE~'>.  xODELT  ^       ^ 

iZZ  AF-    HF1,rO]    01-' e.- ^B*  .^<')^    0)'t-.xr,AtiF 

[3^"  AH*-<AB1.ro1    0)-<a(AF1  ,  fOT    O.M<-.xGA*1B 

r35*  NDIHt-M'xNf^ 

r3'.j  APTO       Tr.  q.NpTM    NE<  ^  ■>?."  >  +  •  " '^NDIH  .  MP  "•  pAF 

r37l  APfO'  ,  V'a(NDiH,NE':»pKPT>  +  -X(NDIri,NB)pAB 

[381  AA«^AC,(;-o.5D    AB 
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FF2DLEV 


I"  "'I 
C133 


'21 


'Ur;i 


.1-,  L  ^     !■■ 
[O'^i    TO 


'-•■fcvr 


CIO-;' 
Ph«-*iVLEV    r 

r-' iK  / 1  1  J 

V    K«-:;!ZE    COVF    F;C:0,K,L,_,k  .Ri 


n    FUrlCTION    ^0    uEinEkmTE    TEKfi; 
ft       _       _ 

F---rtH/^   f 

210* 

L- 


LOO 

-».■!■ 


1 1'" 


■  J  ,  1  M      1-1  > 


rn    jr    rnvARIAMCE    rCN    FOk    IhHuE 


;   K— :• 

2    h,TL,K"l^f  ■^'''--'C  (L    K  )-t-i")  xf  .-(:_    |<  )  :■  4.F  >  •>  -  x  '  .jF 
L,  K  "-'.+/*/•  '.  L,-^y4'F.'X  i,  <-(.L.-K"i  >+F'.  )-x/?F 
>(.  ^K-M  )/L'J0l^'2 

-»■•■  E>L*-L  +  i  ;/lOC'F'» 

Kf-To    0    1    ^.ftl  >  .  R 
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APPENDIX  B 

DERIVATION  OF  THE  POWER  SPECTRUM  AND  AUTOCORRELATION 
FUNCTION  FOR  THE  TWO  POLE  AUTOREGRESSIVE  MODEL 


Power  Spectrum 

|H(e^^)|2   =   HCeJ*^)  •  (H(eJ'^))*   =   H(e>)H(e-3^) 
In  this  case: 


HCeJ'^)   =  ^       H(e-J^)   = 


Calculating  H(e>^  ) H(e"^^) : 

1  11 


(1+ae  ^^) (1+ae^^)  1+ae  ^^+ae^^+a^       l+2acos (w)+a^ 

The   final    result   is: 

SyCco)       =       |H(eJ^)|2      =    ^ (B.l) 

(1+a  )+2acos(oo) 

Autocorrelation  Function 


Starting  with  H(z)  for  this  case: 

H(z)   =  ^-^ 

1  +az 
Per  Ref.  7:p.  158: 

Z~l[H(z)]  =  h(n)  =  (-a)"-u(n)     for  a  <  l 

(u(n)  is  the  unit  step  function) 

Per  Ref.  5:pp.  391-395,  for  the  white  noise  input  case: 

00 

Ry(l)       =  I      h(n)-h(n-£)  (B.2) 

n=-oo 

Substituting  h(n)  above  into  Eq.  B.2 
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(-a)  n=£ 
The  summation  term  may  also  be  expressed  as: 

OO  CO  £-1 

^  (-a) 2^  =  I    (-a) 2^-  ?  (-a) 2^ 
n=£  n=0        n=0 

Per  Ref.l2:p.  8,  the  summation  terms  on  the  right  are  equal 
to: 

y  (_e,)2n  .  1     (,,<i)      I    (-a)2n  .  lj:±:^  (a  <  1) 

n=0  l-(-a)  n=0  1  -  (-a) 

As  a  result: 

V  /  ^2n        1      1  -(-g)    _    (-«)     /  <  1) 

)  (-a)    =  J- T~    ~  2      ^^  ^  ■^' 

n=£  1  -  (-a  )   1  -  (-a)       1  -  (-a) 

Substituting  and  using  (-a)^  =  a^  yields: 

25/  ^ 

RyU)       =     -^--^^^    =     -^^     £lOanda<l  (b.3) 

(-a)        1  -a  1  -a 


105 


APPENDIX  C 

GRAPHICAL  RESULTS  FOR  THE  POWER  SPECTRUM  AND  AUTOCORRELATION 
FUNCTION  FOR  THE  TWO  POLE  AUTOREGRESSIVE  MODEL 


Power  Spectrum 
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3.15 


U3 

/^ 

CS 

/           \ 

^ 

/               \ 

h-  — - 

/                    \ 

•-^ 

v'                                                      ^v 

Z. 

^^^^                                  \^^^ 

O 



< 

o- 

I 

-3.ia 


0.00 

OMEGA 


a.15 


ALPHA=0.5 


-3.16  0.00 

OMEGA 

ALPHA=0.8 


3.16 


-3.16  0.00 

OMEGA 

ALPHA=-0.25 


3.1S 


-3.13 


OMEGA 
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ALPHA=0.75 
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APPENDIX  D 

DERIVATION  OF  THE  POWER  SPECTRUM  AND  AUTOCORRELATION 
FUNCTION  FOR  THE  FOUR  POLE  AUTOREGRESSIVE  MODEL 


Power  Spectrum 

|H(eJ^)|2      =      H(e>)-H(e">) 
In  this   case: 
H(e">)       = 


H(eJ'^)       = 


(l-ae  -*  Q-'   )  (1-ae-^  e     )  1-ae-^  -ae  -'  -fa  e-^ 

1  1 


11        ~J6  ~JWv  ,n        j6  -j(jo>  ,        -j  (9-Ko)        j  (6-co)  ,    2  -j2co 

(l-ae  -'  e  -^   )  (l-ae-'  e  -"   )  l-ae  -^  -ae-'  -i-a  e  -' 

Multiplying  the  above  expressions  yields: 
H(e>)H(e-J")   = 


j(e-co)      -i(e+to)    2  -j2co     i(e+oj)    2  j2e    2 

l-ae-"  -ae  -^  -fa  e  -^     -ae-^  -fa  e-"     +a 

3  j(0-co)        -j(e-w),    2,    2  -120     3  -jO-co) 
-  a  e-^  -ae  -'  +a  +a  e  -*     -a  e  -^ 

.      2  j2oj     3  j  (0+co)      3  -j  (0+co)  ,    4 
+  a  e-^     -a  e-^  -a  e  -^  +a 


Combining  terms: 
H(e3W)  .H(e-J^)   = 


T^^   2^  4     3,  -j(0-aj)_^j(e-co),^  2,   j20^-j20. 
l+2a  +a  -a    (e  -•  -fe-"  )  +a    (e-^     +e  -^      ) 

3,   j  (0-Hjj)  ,    -j(0-fco)v,    2,   j2go,    -j2a), 

-  a    (e-^  '  -fe  -^ '         )  +a   (e-'     -fe  -'     ) 

-   j(0-aj),    -j(0-oo)>       ,   j  (0+aj)^-j  (0+co). 

-  a  (e-"  +e  -^  '         )  -a  (e-*  '  +e  -'  ) 


Using  Euler's  relation  and  combining  terms: 


HCeJ*^) -HCe'J^)       = 


4         2  3  2 

1+a  +2  (a  -cos(0-oj)  [a  +a]+a  cos  (20)-cos  (0+aj) 


3  2 

X    [a  +a]+a  cos(2to) ) 
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Using  cos  (6-0))  +  cos(9+co)    =    2cos  (  6  )  cos  (  co  )  ,   and  since 
Sy(co)  =  H(e3^)  •H(e":i'^)  •  o2,  with  o^    =   1   the  final  result  is: 


Sy(w)       =     T 5 ^ 2 (^-1) 

1+a  +2  (a  -2  [a  -a] cos  (9)  cos  (co)  +a    (cos  (26)  +cos  (2aj) ) ) 


Autocorrelation   Function 


H(z)       = 


-,        je  -1       -je  -1.    2-2  ,     ,   j6^  -j6,    -1^   2  -2 

1-ae-^   z     -ae  -^   z     +a  z  1-  (e-^  -fe  -'    )  z     +a  z 

Using   Euler's   relation: 


H(z)       =    zi — 2~^  1^1    ^    ^^ 

l-2acos(e) z     +a  z 


Per  Ref.  7:pp.  204-216,  partial  fraction  expansion  can  be 
used  to  find  the  inverse  Z  transform.  To  do  so  H(z)  can  be 
expressed  in  the  form: 

H(z)   = 


,.      -je  -i>  ,,   je  -1,    ,   -j6>  -   je, 

(l-ae  -^  z  )  (l-ae-'  z  )     (z-ae  -^  )  (z-ae-*  ) 

Using  the  partial  fraction  expansion  and  table  look  up  [Ref 
7:p.  158]  yields: 


n 

h(n)   =   ■  ^,^.  cos(ne+9-^)  •u(n)   for  a<l  (D.2) 

sin  (6)         2 


Since  cos(e  -  y)  =  sin(e),  the  final  expression  for  h(n)  is 


n 

h(n)   =   "  ■^,  sin((n+l)e)  •u(n)      a<l  (D.3) 

sin(e) 
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For  simplicity  in  further  derivation  of  Ry(ji)  based  on  h(n), 
Eq.  D.2  will  be  used. 

Using  the  expression  for  the  autocorrelation  function  of 
a  random  process  represented  by  the  above  filter  with  a 
white  noise  input  [Ref.  4:pp.  391-395]: 


00  -£       oo 

Ry(^)  =    I     h(n).h(n-e)  =  —^ I     a  •oos(ne+e-^) 

n=-~  sin  (6)  n=£ 

•  cos((n-£)e+e-^)  (D.4) 


n  =  5,  in  the  summation  index  since  h(n)  is  causal.  £  is 
assumed  to  be  greater  than  zero  here.  For  il  <  0,  Ry(£)  = 
Ry(-£)  by  symmetry  of  the  autocorrelation  function  [Ref. 
5:p.  388],  so  we  can  proceed  assuming  only  positive  values 
of  £. 

Using   the   trigonometric   identity   for   a   product   of 
cosines: 


-£         °° 

Ry(£)       =     —^ I     a^^i  .cos(£e)+i^s(2ne-£e+2e-7T)]      £>0 

sin^(O)   n=£ 


—  £  oo  -£  oo 

R   (^)      =     cos{iQ)'a       .      ^     „2n^_a    ^     a^^.  (2ne-£  0+26-7:)         (D.5) 

^  2sin   (6)  n=£  2sin   (9)   n=£ 

£  >  0 

CO  °o  £-1 

Using  I        ^         1      ~      1         ^^^   standard   geometric  progression 

n=£      n=0    n=0 
identities  [Ref.  ]2:p.  8]: 

£-1         ,  ^£ 

)   ct    =   T and      )   a    =   -,  _^ 

nio         ^""  n=0         ^  ^ 
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For  the   first  term   in  Ryrd) 


2n  V       2n       ''v'     2n     _         1  l-a^^ 


Z  a         =       I     a       -     I 


a 


~  2  2 

n=£  n=0  n=0  1-a         1-a 


2£ 

a 


l-a^ 


(D.6) 


For  the  second  term  in  Ry(jl): 

let        4>  =  -  5,6  +  2Q  -  7T 

Using  Euler's  relation: 

j  (2ne+(t))  ^  -j  (2ne+({)) 


cos(2ne+(t))  = 


V   2n    ,o  n^A^     1  V   2nr  i(2ne+(j))  ,   -j(2ne+(})), 
I     a     •cos(2ne+(j))  =     j    I     a      [e-'  ^+e-'     ^J 

n=il  n=£ 

^  n=il  ^  n=il 

For  large  n,  it  is  evident  that  the  a^n  term  will  tend 
to  make  the  term  in  each  sum  approach  0  for  a  <  1 ,  and  thus 
ensures  convergence  and  a  closed  form  expression  for  each 
term.  Pursuing  the  mathematics  required  to  find  this  closed 
form  expression  we  have: 

2  n=£  n=i  n=0         n=0 
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1  r       2n     j  (2ne+(|))      ^     e^^  1  _  l-Cae^^)^^, 

2  ^0  2    %     ,     je,2       ,     .     je,2   ^ 

2£     j  (2eJl+(t)) 


2(l-(ae^^)2) 


(D.8) 


For   the   conjugate   term  we  must   have: 


1      °°     o,.        •  /T   Q.A^  2£     -j(2e£+(})) 

n=Jl  2  (l-(ae  -^   )    ) 


Next  a  common  denominator  must  be  found  to  sum  these  two 
terms: 


21    j(2e£+4))   ,,  .  -je.2. 

i     .e-"     ^    (l-(ae  -*  )  ) 


2(l-(ae^®)^)    (l-(ae  ^^)^) 


21   j(2e£+(}))_  (2£+2)  j  (2e£-2e+(})) 
2  (1-  (ae^  )  -  (ae  -^  )  -i-a  ) 


a 


21     -j(2e£+(|))    ,T  ,  j0,2, 
•e  -'     ^'        (l-(ae-'  )  ) 


2(l-(ae"^®)^)    (l-(ae^®)^) 


2£  -j(2eil-k}))_  (2£+2)  -j  (2e£-2e+(i)) 

2(l-(ae-'   )    -(ae  -"   )   +a  ) 
By  Euler's   relation: 

l-Cae-"   )   -(ae  -"   )   +a       =     l-2a  cos(29)-i-a 

Adding  the  terms  with  the  common  denominator  yields 

2£  j(2e£+(j))  ^2  i(2e£-20+(}))_^-j(2e£+(t))  ^2  -j(2e£-2e+(^), 
a   [c        —a  e  "rs         —a  s  j 

2(l-2a^cos(20)+a^) 
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Again  using  Euler's  relation  the  above  expression  reduces 
to: 

g^^  [cos  (29£+(J))  -a^cos  (26£-26+(|))  ]  .^   .^. 

2        4  \u.iU) 

l-2a  cos(29)+a 

which  is  the  sum  of  the  last  two  terms  in  Eq.  (D.7). 

Substituting  Eq.  (D.IO)  and  Eq.  (D.6)  into  Eq.  (D.5)  yields: 


R   iZ)     =     cos(£9)-a"^^a^j 
^  2sin    (6)        1-a 

-£  2i  2 

a  rCi      [cos(295,+(}))-a  cos(29£-29+({))  J  /n    ^1^ 

+  2 L 2 4 ^LJ.i±; 

2sin   (9)  l-2a  cos(29)+a 

Combining  and  canceling  terms  and  substituting  for  <t>   and 
noticing  that  the  same  result  must  hold  for  £  <  0  we  have: 


I  11  II  2  I     1 

,».      _  g        fOOs{\i\Q)    ,  cos(  (24-|  £[  )9-7T)-a  cos(  I  £[  9-7t) -■  /^    -12^ 

^  2sin    (9)      1-g  l-2g  oos(29)+g 
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APPENDIX  E 

GRAPHICAL  RESULTS  FOR  THE  POWER  SPECTRUM  AND  AUTOCORRELATION 
FUNCTION  FOR  THE  FOUR  POLE  AUTOREGRESSIVE  MODEL 
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APPENDIX  F 

DERIVATION  OF  THE  POWER  SPECTRUM  AND  AUTOCORRELATION 

FUNCTION  FOR  THE  FOUR  POLE  AUTOREGRESSIVE  MODEL 

(WITH  TWO  POLES  ON  THE  REAL  AXIS^ 


Power  Spectrum 

Using  initial  results  from  Appendix  D  with  the  necessary 
modifications  (including  9  =  0)  we  have: 


H(e-:^)       = 


„-'       ^r\]      z^-'       -uru      r\i       d-' 

a' 


1IUJ  JUU 

-a.e-^   ~(\Q     +ct  a  e 


HC®-"^)       -      ,  -jco         -jw  -j2ca 

1-a  e  -"  -C6.  e  -"  +a  a,  e  -* 


Therefore: 


H(eJ^)-H(e-J^)      = _.^  -^ ::^ ja.^2^ 

1-a  e       -a.  d  -'  +a  a,  e  -'     -a  e-^  ""^a     a°'b 

2       -iw         io),  ,2         2  -joo 

i2co     2        jo)  2  ju),    2  2 

+  a  ex.  e-'     -a  cl  e-"  -a  cCe-'  ^ctgO^ 

Combining  terms  and  using  Euler's  relation: 


H(eJ^)-H(e-J^)   =   ^ ^ 2 

l-(a  H-a+a  a+a  a  )  •2oos(w)+2a^o(j^cos(2(ja)-i-a^ 


2  2  2 


4-  2a^(^+(^+a^a^ 
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Assuming  a^  =  1   and   since   Sy(a3)  =  H(el)^) -HCe"!)^  )  •  o2 ,   the 


final  result  is: 


Sy(co)   —  2     2  2 

1-2  (a  +a +a  a, +a  a,  )cos(a3)+2a  a  cx)s(2aj)+a  +2a  a, 

+  CL +a  a,  ( F .  1 ) 

Autocorrelation  Function 


H(z)   = 


,   je  -1   -j0  -1    -2 

l-a  e-^  z   -a,  e  -^  z  -i-a  a,  z 
a       D        a  D 


Letting  6=0 


1  2 

H(z)   =   1  ^ 


-I      -2  ~   2 

1-  (a  +a,  )  z  -i-a  a,  z       z  -  (a  +(V )  z+a  a. 


z  >  a 


Expanding  in  terms  of  partial  fractions  we  have; 


a.      _     a^ 


H(z)  =  -^—  '  ~^~  + 


V^   ^-^'a   °t-^a   ^-"b 


This  corresponds  to  the  impulse  response 


a 


n   .      '\>  n. 


h(n)   =   ( -a  +  •  a,  )  -uCn)    a  <  1  fF  2^ 

a  -a,    a   a -a  lo'  K^-^) 

a  D       Da 

Proceeding  as  in  Appendix  D: 
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n+1  n+1       n+l-J.       n+l-l 


^  n=-oo  nU  V^       V«a       V%     V 


•a 
a 


{I  >  0) 


2n+2-i  n+1     n+l-£  n+l-£     n+1 

a 
a  ,     a 


°°     a  a       .a,  a  •a. 


n=£    (a^-Oj^) 


2n+2-£ 
+  — 2^  (£>0)  (F.3) 


(c^-a^) 
Since      (ct  -a  )   =  -(a -a  ): 


K-%^^     =      ^V^a^^     =     -K-%^^%-%^     =     -(V^'a^^V^^ 


So: 


„    ,„,  1        r  2-£     V       2n  l-£     V    /  nI^ 

R  (£)  =    T-[a  /  ct  -  a  a  /    (a  a,  ) 

y  ,  \2     a         ^n     a  a  D         -^^     a  d 

^  (a  -a  )  n=ic  n=i 


00 


l-£     V    /         \r^   ,     2-£     V     ^n,        ,„   ^  n\ 

Va         V^"a)  ^         ^a"^^         ^     "^ 

n=l  n=£ 


Continuing  with  the  same  principles  and  assumptions  as  in 
Appendix  D,  we  have: 


V   2n     V   2n    V   2n      1    1-ct      a 
la         =       la       -la 


~     2       2       2 

n=£        n=0      n=0        1-a     1-a     1-a 


£-1         ,  T   £  £ 

V  n      r   n    V  ri      1  1-a  a 

^.                  ^f.               ^f.                 1-a  1-a  1-a 
n=£      n=0     n=0 
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Making  the  appropriate  substitutions  in  the  expression  for 
Ry(il)  ,  we  have: 

%''^^  ,  ,2^°'a  \      T       ^a"b  ^1-a  a  ^ 

•^       (a  -a.  )      1-a  a  d 

A/     A/  ^  A/ 

a  D        -'-"^ 

Combining  terms  yields  the  final  expression: 

2+1         £+1      £+1    2+£ 
,    a     a   a, -i-a  ex,     a, 
_  ,  .         1   r  ^ a  ^  a  D    _b , 

^y^^^   "   ,     ^2^,   2      1-a  a       i   2  ^  (F.5) 

(a  -a,  )  1-a        a  d      -^"^ 
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APPENDIX  G 

GRAPHICAL  RESULTS  FOR  THE  POWER  SPECTRUM  AND  AUTO- 
CORRELATION FUNCTION  FOR  THE  FOUR  POLE  AUTO- 
REGRESSIVE  MODEL  rWITH  TWO  POLES  ON  THE  REAL  AXIS) 

Power  Spectrum 
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APPENDIX    H 
LAPLACIAN    INVERSE    FILTER    FORM.q 
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APPENDIX    I 

CONVOLUTION   OF    LAPALACIAN    DIFFERENCE    OPERATOR   AND 
VARIOUS    SIZE    FIR    INVERSE    FILTERS 
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